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We derive relations between various observables for N particles with zero-range or short-range 
interactions, in continuous space or on a lattice, in two or three dimensions, in an arbitrary external 
potential. Some of our results generalize known relations between large-momentum behavior of 
the momentum distribution, short-distance behavior of the pair correlation function and of the 
one-body density matrix, derivative of the energy with respect to the scattering length or to time, 
and the norm of the regular part of the wavefunction; in the case of finite-range interactions, the 
interaction energy is also related to dE /da. The expression relating the energy to a functional of the 
momentum distribution is also generalized, and is found to break down for Efimov states with zero- 
range interactions, due to a subleading oscillating tail in the momentum distribution. We also obtain 
new expressions for the derivative of the energy of a universal state with respect to the effective 
range, the derivative of the energy of an efimovian state with respect to the three-body parameter, 
and the second order derivative of the energy with respect to the inverse (or the logarithm in the 
two-dimensional case) of the scattering length. The latter is negative at fixed entropy. We use 
exact relations to compute corrections to exactly solvable three-body problems and find agreement 
with available numerics. For the unitary gas, we compare exact relations to existing fixed-node 
Monte-Carlo data, and we test, with existing Quantum Monte Carlo results on different finite range 
models, our prediction that the leading deviation of the critical temperature from its zero range 
value is linear in the interaction effective range r e with a model independent numerical coefficient. 

PACS numbers: 



I. INTRODUCTION 



The experimental breakthroughs of 1995 having led to the first realization of a Bose-Einstein condensate in an 
atomic vapor [IHH have opened the era of experimental studies of ultracold gases with non-negligible or even strong 
interactions, in dimension lower or equal to three I n these systems, the thermal de Broglie wavelength and the 

mean distance between atoms are much larger than the range of the interaction potential. This so-called zero-range 
limit has interesting universal properties: Several quantities such as the thermodynamic functions of the gas depend 
on the interaction potential only through the scattering length a, a length characterizing the low-energy scattering 
amplitude of two atoms. 

This universality property holds for the weakly repulsive Bose gas in three dimensions [9J up to the order of expansion 
in (rta 3 ) 1 / 2 corresponding to Bogoliubov theory jlOj . n being the gas density. It is also true for the weakly repulsive 
Bose gas in two dimensions [llT - ll3| . even at the next order beyond Bogoliubov theory For a much larger than 
the range of the interaction potential, the ground state of N bosons in two dimensions is a universal iV-body bound 
state [15l - ll9j |. In one dimension, the universality holds for any scattering length, and the Bose gas with zero-r ang e 
interaction is exactly solvable by the Bethe ansatz both in the repulsive case [20| and in the attractive case [2l|, |22[ . 

For spin 1 /2 fermions, the universality properties are expected to be even stronger. The weakly interacting regimes 
in 3D 12314271 and in 2D [28| are universal, and the ID case is also solvable by Bethe ansatz for an arbitrary interaction 
strength [29[|30j]. Universality is expected to hold for an arbitrary scattering length even in 3D (see however l3ll), a s 
was recently tested by experimental studies on the BEC-BCS crossover using a Feshbach resonance, see e. g. [a. |32h39( | . 
and in agreement with unbiased Quantum Monte Carlo calculations [4Ch45| | . A similar universal crossover from BEC 
to BCS is expected in 2D when the parameter ln(fci?a) varies from — oo to +00 (4fj| - [5ci| . Universality is also expected 
for mixtures in 2D 50-5j|, and in 3D for Fermi- Fermi mixtures below a critical mass ratio [5ll [53l l54j|. 

In the zero-range regime, it is intuitive that the short-range or high-momenta properties of the gas are dominated 
by two-body physics. For example the pair distribution function g( 2 )(ri 2 ) of particles at distances r i2 much smaller 
than the de Broglie wavelength is expected to be proportional to the modulus squared of the zero-energy two-body 
scattering wavefunction <^(ri2), with a proportionality factor A g depending on the many-body state of the gas. 
Similarly the large momentum tail of the momentum distribution n(k), at wavevectors much larger than the inverse 
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de Broglie wavelength, is expected to be proportional to the modulus squared of the Fourier component of the zero 
energy scattering state 0(k), with a proportionality factor A„ depending on the many-body state of the gas: Whereas 
two colliding atoms in the gas have a center of mass wavevector of the order of the inverse de Broglie wavelength, 
their relative wavevector can access much larger values, up to the inverse of the interaction range, simply because the 
interaction potential has a width in the space of relative momenta of the order of the inverse of its range in real space. 

For these intuitive reasons, and with the notable exception of one-dimensional systems, one expects that the mean 
interaction energy E- mt of the gas, being sensitive to the shape of g^ at distances of the order of the interaction range, 
is not universal, but diverges in the zero-range limit; one also expects that, apart from the ID case, the mean kinetic 
energy being dominated by the large-momentum tail of the momentum distribution, is not universal and diverges 
in the zero-range limit, a well known fact in the context of Bogoliubov theory for Bose gases and of BCS theory for 
Fermi gases. Since the total energy of the gas is universal, and E- lnt is proportional to A g while -Ekin is proportional 
to A n , one expects that there exists a simple relation between A g and A n . 

The precise link between the pair distribution function, the tail of the momentum distribution and the energy of the 
gas was first established for one-dimensional systems. In [20] the value of the pair distribution function for m = was 
expressed in terms of the derivative of the gas energy with respect to the one-dimensional scattering length, thanks 
to the Hellmann-Feynman theorem. In [55| the large momentum tail of n(k) was also related to this derivative of the 
energy by using a simple and general property of the Fourier transform of a function having discontinuous derivatives 
in isolated points. 

In three dimensions, results in these directions were first obtained for weakly interacting gases. For the weakly 
interacting Bose gas, Bogoliubov theory contains the expected properties, in particular on the short distance behavior 
of the pair distribution function [56l458| | and the fact that the momentum distribution has a slowly decreasing tail. 
For the weakly interacting two-component Fermi gas, it was shown that the BCS anomalous average (or pairing field) 
(^-|-(ri)^(r 2 )) behaves at short distances as the zero-energy two-body scattering wavefunction 0(ri 2 ) b 59j, resulting 
in a function indeed proportional to |</>( r i2)| 2 at short distances. It was however understood later that the 
corresponding proportionality factor A g predicted by BCS theory is incorrect [o0| . e.g. at zero temperature the BCS 
prediction drops exponentially with 1/a in the non-interacting limit a — > CP, whereas the correct result drops as a 
power law in a. 

More recently, in a series of two articles (HH, [62]] , explicit expressions for the proportionality factors A g and A„ were 
obtained in terms of the derivative of the gas energy with respect to the inverse scattering length, for a two-component 
interacting Fermi gas in three dimensions, for an arbitrary value of the scattering length, that is, not restricting to 
the weakly interacting limit. Later on, these results were rederived in [63l - l65j . and also in [66[ with very elementary 
methods building on the intuition that g^ oc |</>(Vi2)| 2 at short distances and n(k) cx |</>(fc)| 2 at large momenta. 
These relations were recently tested by numerical four-body calculations [67| . An explicit relation between A g and 
the interaction energy was derived in 65 1. Another fundamental relation discovered in [6l[ and recently generalized 
in [68[ to bosons, to Fermi-Bose mixtures and to fermions in 2D, expresses the total energy as a functional of the 
momentum distribution and the spatial density. 

In the present work we derive generalizations of the relations of [2(| [H, [U [62|, [fjjj [68[ to two dimensional gases, 
to the general case of a mixture of an arbitrary number of atomic species and spin component, and to the case of a 
small but non-zero interaction range (both on a lattice and in continuous space). We also find entirely new results for 
the first order derivative of the energy with respect to the effective range and, in presence of the Efimov effect, with 
respect to the three-body parameter, as well as the second order derivative with respect to the scattering length. 

The article is organized as follows. In Section|TT]we treat in detail the case of two-component Fermi gases. Relations 
holding for any system eigenstate for zero-range interactions are derived in Section III Bl and summarized in Table [II] 
We then consider lattice models (Tab. Hill Sec. Ill C|) and finite-range models in continuous space fTab. HVl Sec. Ill Dl) . 
In Section III El we derive a model-independent expression for the correction to the energy due to a finite range or a 
finite effective range of the interaction. The generalization to thermodynamic equilibrium, where the system is in a 
statistical mixture of eigenstates, is discussed in Section lll Fl In Section HTT1 we turn to the case of spinless bosons. We 
focus on the case of zero-range interactions where, in 3D, the Efimov effect leads to modifications or even breakdown 
of some relations, and to the appearance of a new relation. Then we show briefly in Section HVl how to treat the case 
of an arbitrary mixture and present results for zero-range interactions (Tab. IVI[) . Finally we present applications of 
exact relations: For three particles we compute corrections to exactly solvable cases and compare them to numerics 
(Sec. IV Ap . and we check that exact relations are satisfied by existing fixed-node Monte-Carlo data for correlation 
functions of the unitary gas. We expect from our expression for the leading finite-range correction to the energy that 
the leading finite-range correction to the critical temperature in the BEC-BCS crossover depends only on the effective 
range of the interaction, an expectation that we test against the Quantum Monte Carlo calculations of [4l|, Ejj . We 
conclude in Section IVT1 
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II. TWO-COMPONENT FERMIONS 

In this Section we consider spin-1/2 fermions. For a fixed number N„ of particles in each spin state a =t,i, one 
can consider that particles 1, . . . , jV+ have a spin f and particles iV+ + 1, . . . , jV+ + N\ = N have a spin J,, i.e . th e 
wavefunction "0( r ij • ■ • , r _/v) changes sign when one exchanges the positions of two particles having the same spin [l32j . 

A. Models 

Here we introduce the three models used in this work to model interparticle interactions. 

1. Zero-range model 

In this well-known model (see e.g. [69MT6| and refs. therein) the interaction potential is replaced by contact 
conditions on the many-body wavefunction: For any pair of particles i ^ j, there exists a function A i j 1 hereafter 
called regular part of ip, such that in 3D 

t/)(ri,...,r N ) =^(—--)Ai j (n ij ,(T k ) k jt i j) + 0(ri j ), (1) 



and in 2D 



a 



V»(ri, . . . ,rjv) = ln(nj/a)Aij (Rij^r^kjti j) +0(nj), (2) 

ru-H) 



where the limit of vanishing distance r,-j between particles i and j is taken for a fixed position of their center of mass 
Rjj = (rj + Yj)/2 and fixed positions of the remaining particles {Yk)k^i,j- Fermionic symmetry of course imposes 
Aij = if particles i and j have the same spin. When none of the r^'s coincide, there is no interaction potential and 
Schrodinger's equation reads 

H1>(r 1 ,...,T N )=Ei/)(T 1 ,...,T N ) (3) 

with 



N r t2 



-^ Ti+ U{v i ) 



tP (4) 



where m is the atomic mass and U is an external potential. The crucial difference between the Hamiltonian H and 
the non- interacting Hamiltonian is the boundary condition (|H2p . 



2. Lattice models 

These models were used for q uantum Monte-Carlo calculations pol - l43l . liH W?\. They can also be convenient for 
analytics, as used in 0, [78|, l79j | and in this work. Here particles live on a lattice, i. e. the coordinates are integer 
multiples of the lattice spacing b. The Hamiltonian reads 



H = H a +g W (5) 



where 



N r .2 



-|-A rs +C/(r 4 ) 
2m 



(6) 



W = Y, S ^b- d (7) 
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Three dimensions 


Two dimensions 


ip(n,...,r N ) = [—--) A lJ (R i:j ,(r k ) k7 , l j) + 0(r tJ ) 


i/>(ri, ...,r N ) = \a{rij/a) A tj (Ry, (r k )k&,j) + 0{r l3 ) 





TABLE I: Notation for the regular part A of the many-body wavefunction appearing in the contact conditions (first line) and 
for the scalar product between such regular parts (second line). 



Three dimensions 



dE Aufi 



d(—l/a) m 



(A, A) 



n— v 7 4 n ^ 4yrm dE 
C = hm k n<j (k) = — — 

fc^+oo h 2 d(— la) 



I 



+ 2' 2^ r^o (47r) 2 r 2 



-E — -Etr. 



fi 2 C 
47rma 



+E 



r d\ 


ft 2 fc 2 


n CT (k) - 


C" 


1 (2tt)3 


2m 


" fc 4 _ 



/ 



i — 1 cr 

C m 
~4^ r ~ 3h? 



2 2 / r^O 

ri 2 c* \ 



TV 



E - Strap - r + o(r 2 ) 



1 d 2 E ( 4irh 



2d(-l/a) 2 V m 



n,E n j£E 



E — E n 



d 2 F \ f d 2 E \ 

d(-l/a) 2 ) T <0 > U(-l/a)V s < 



^ = jg£d(-l/a) /A 
dt 47rm dt V ; 



Two dimensions 



d£ 2tt/1 2 



d(lna) m 



(A, A) 



4 /, v 2™ dE 

C = hm k n a [k) = —pr~n\ \ 

fc->+oo h 2 d(lna) 



d 2 Rg^ (R+|,R- 



2/ r^o (2tt) 



C . 2 

In r 



E — Strap = lim 



+E 



h 2 C , /aAe 7 

m 

27rm V 2 

d 2 fc /i 2 fc 2 x _ 
r(k) 



fc<A (2^) 2 2m 



2 fl<&> (r+£,R-£) = N, + f- r 2 In r + 0(r 



2 / r^O 



TV 



C 
--. — ' 

47T 



(£- Strap) r 2 + o(r 2 " 



2h 2 



I d 2 S _ (2-Kh 2 \ 2 y, (,4 (n) , A)| s 



2d(lna) 2 V ™ 7 ^ £ - £ 

n,E„^E 



d 2 F \ ( d 2 E , 

< 0, w < 



d(lna) 2 



Vd(lna) 2 



dE 

"<t7 



fi 2 C d(lna) / . . 



TABLE II: Relations for two-component fermions with zero-range interactions. The regular part A is defined in Tabled] Lines 
1-7 hold for any eigenstate, and can be generalized to finite temperature by taking a thermal average in the canonical ensemble 
and by taking the derivatives of E with respect to a at constant entropy S. Line 8 holds in the canonical ensemble. Line 9 
holds for any time-dependence of scattering length and trapping potential and any corresponding time-dependent statistical 
mixture . 
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G ■: 



4nm dE 
h 2 d(-l/a) 



C ■ 



2nm dE 
H 2 d(lna) 



dE 47r/T 



d(— 1/a) m 



(A, A) 



dE 



d(ln a) m 



(A, A) 



E iB 



h 2 \ 2 C 



ml g 



E — -Etrap — 



h 2 c 

Anma 



E 



d 3 k 

(27T) 



n CT (k) - C* 



/ ft 2 



E - -Etrap = lim 



h 2 C ( aqe 



■In 



+E 



d 2 k 



(2tt; 



7'k 



o 2-Km 



n a (k)-C-^—T- t lr 



2mek 2m(ek — e q ) 



u^l« \{A in) ,A)\ 2 

n,E nl tE 



d 2 F \ <Q /d 2 ^ () 



E 



6 3 5 | 2) (R,R) = ^|0(O)| 2 



E & Mi( R ' R )=(^l^°) 



In the zero-range regime k typ b <C 1 



^ b*$ (R + £,R - £) * ^|<^(r)| 2 for r « fc^ 



(4tt) 



£ & 2 5 } 2 >(r+I,r-I 



^|0(r)| 2 forr«fc- 1 p 



/ ft 2 \ 2 
n CT (k) ~ C ( J for k » fc typ 



TABLE III: Relations for two-component fermions in a lattice model. C is denned in line 1. 



in first quantization, i.e. 

H « = £ / 4%e k 4(k)c,(k)+^6 d C/(r)(^^)(r) (8) 
= Y, bd (^A^)(r) (9) 

r 

in second quantization. Here d is the space dimension, e k is the dispersion relation and c£.(k) is creates a particle 
in the plane wave state |k) defined by (r|k) = e ikr for any k belonging to the first Brillouin zone D = (— f , f ] d - 
Accordingly the operator A in is the discrete representation of the Laplacian defined by — ^(r|A r |k) = ek(r|k). 

The simplest choice for the dispersion relation is ek = — — [2 EH EE IZSi lzl| • Another choice, used in [4l|, [z3 ; is 
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Three dimensions 



Two dimensions 



C : 



4nm dE 
h 2 d(-l/a) 



C : 



2irm dE 
H? d(lna) 



E- m 



C 



(4tt) 



d 3 rV(r)\4>{r)\ 



Eint = j^r I f/L, '' r( '' , i" ( '' , i L> 



E 



-E — i?tr 



h 2 c 

4-Kma 



-B - Strap = lim < In — 

i?->cx) liama \ a 



(2tt) 3 2m 



C 



nCT(k) ~ (4^ l0 ' (fc)l 



E 



(2tt) 2 2m 



n CT (k)-^|^(fc)| 



In the zero-range regime k typ b <C 1 



'Rgfl ( R+ £,R- £) ^ ^|^>(r)| 2 for r « Ky \ 



(Air) 



2 Rg^ ( R+ £,R- £) c ^Wr)| 2 for r « fc^ 



(27T) 



C ~ 2 



C 

lo-(k) ~ 7^^l^( k )| 2 forfc>fcty P 



TABLE IV: Relations for two-component fermions with a finite-range interaction potential V(r) in continuous space. C is 
defined in line 1. 



h 2 

the dispersion relation of the Hubbard model: ek = — 7 [1 — cos( hb)]. More generally, what follows applies to 

too ' 

i=l 

ft 2 fc 2 

any ek such that ek — > sufficiently rapidly and e_k = £k- 

6— s-o 2m 

A key quantity is the zero-energy scattering state <j>(r), defined by the two-body Schrodinger equation (with the 
center of mass at rest) 

^ A r +5o^W)=° ( 10 ) 



rn 



and by the normalization conditions 



<j)(r) ~ --- in 3D (11) 

r>b r a 

<p(r) ~ ln(r/a) in 2L>. (12) 

A straightforward two-body analysis, detailed in App. [3] yields the relation between the scattering length and the 
bare coupling constant go- 

1 to f d 3 k 1 „ . . 

in 3D (13) 

g 9^0 27m 2 (2tt) 2 2(e q - e k ) 

where 7 = 0.577216 ... is Euler's constant and V is the principal value. Other useful relations derived in App. lAl arc 

, . Airti 2 , . 

0(0) = in 3D (15) 

mgo 

0(0) = in2L> (16) 

mg 



go 47r?i 2 a (27r) 3 2e k 

1 
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Three dimensions 



Two dimensions 



dE 



4-irtT 



d (- 1 / a )J Rt m 



(A, A) 



dE 2nh 2 



d(lna) m 



(A, A) 



C = lim k 4 n(k) = ( Tu~TT \ ' 

k-++oo h? \d(—l/a) J R 



C= lim fc 4 n(fc) = f E 



d 3 Rg (2) [H + ^,R- ^ 



C 1 



2/ r^o (47r) 2 r 2 



>o (2?r) 



O , 2 

m r 



-E — E t r 



h 2 c 



8ivma 
f d 3 k h 2 k 2 r 
+ h m / 7 — 



-E — Etr 



lim 



/i 2 C , / aAe 7 
■ In 



47rm 



d 2 fc fi 2 fc 2 



9 2 E 



2 V<9M/a) 2 



2\ 2 



E 

n,E n ^E 



\{A^ n \A)\ 2 



E — E n 



I d 2 £ _f2irh 2 \ 2 l(^ (n) ,^)| 5 



2d(lna) 2 V m J ^ E - E, 

K ' V 7 n,E n ^E 



x / dC / dr4 . . . drjv |_B(C, r4, . . . , rjv) 



/ d 2 F 
Vd(lna) 



d 2 g 
d(lna) 2 



< 



< 



TABLE V: Main results for spinless bosons in the limit of a zero range interaction. In three dimensions, the derivatives are taken 
for a fixed three-body parameter R t . As discussed in the text, in three dimensions, the relation between energy and momentum 
distribution is valid if the large cut-off limit A — > +oo exists, which is not the case for Efimovian states (i.e. eigenstates whose 
energy depends on Rt). In the last relation in three-dimensions, B is the three-body regular part defined in (|152|l . 



and 

, . . AttH 2 d(—l/a) „ , , 

|0(O)| 2 = L^J in3jD (17) 

|0(O)| 2 = in 2D. (18) 

m dg 

In the zero-range limit (b — > with go adjusted in such a way that a remains constant), the spectrum of the 
lattice model is expected to converge to the one of the zero-range model [Il|,[79j], and any eigenfunction ?p(ri, . . . ,rjv) 
of the lattice model tends to the corresponding eigenfunction of the zero-range model, provided all interparticle 
distances remain much larger than b. Let us denote by l/kt yp the typical length-scale on which the zero-range model's 
wavefunction varies: e.g. for the lowest eigenstates, it is on the order of the mean interparticle distance, or on the 
order of a in the regime where a is small and positive and dimers are formed. The zero-range limit is then reached if 

hypb < 1. 

For lattice models, it will prove convenient to define the regular part A by 



V>(ri,...,r< = R^,...,!^ = Ry,...,rjv) = 0(0) A^Ry, (r fe ) fe ^). 



(19) 
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Three dimensions 



Two dimensions 



dE 2nh 2 



d(-l/«w) |U CT 



{A,A) C 



dE wh 2 



d(lna CT<T ,) fj, a 



{A,A) a 



C a = lim fcV(k) = Y(l + 8 a 



/I 2 fl(-l/a„0 



C CT = lim fc 4 n CT (k) = V(1 + <U0 



A-Kji aa , dE 
h 2 c^lna^, 



-r,R 



as i 



m CT + m 



-r,R- 



m CT + m„i J r->o 



(1 + <U0 



'2-kK 2 d(-l/a aa >) r 2 



X 



i 2 

m r 



7r/l 2 a(lna„ 



E — Etm 



+E 



_L 6 

4-* a„„, d(-l 



dE 

9(-l/a CT0 ./) 



— -Etrap = lim 



E dE ( a aa ,Ke 1 

£ ,9(ln<w) n 



r d\ 


h 2 k 2 


n CT (k) - 


C CT " 


J (2tt)3 


2m a 


" fc4 . 



(2tt) 2 2m^ 



■w CT (k) 



1 <9 2 £„ /27rfi : 



2 0(-l/<w)a 



1 <9 2 £ n / vrfi 



* 2 x 2 |(A ( "' ) ,A ( " ) ) 
2<9(lnaw) 2 V/W7 n , tl f^ En E n -E n , 



d 2 F 



9(-l/<w) 2 



< 



9 2 F 



9(lno, ^ 2 



aa' ) J T 



< 



d 2 E 



9(-l/<w) 2 



< 



a 2 £ 



<9(lna CTCT /) 2 



< 



TABLE VI: Main results for an arbitrary mixture with zero-range interactions. In three dimensions, if the Efimov effect occurs, 
the derivatives must be taken for fixed three-body parameter(s) and the expression for E in line 4 breaks down. 



In the zero-range regime k typ b <C 1, we expect that when the distance between two particles of opposite spin is 
<C 1/fctyp while all the other interparticle distances are much larger than b and than , the many-body wavefunction 
is proportional to 0(r;j), with a proportionality constant given by (|19p : 

■0(ri,...,rjv) - 4>( r j ~ r i) (rk)kjH,j) (20) 

where Ry = (r^ + r 3 )/2. If moreover ^> b, <fi can be replaced by its asymptotic form (|11|12|) : since the contact 
conditions ([T]), ([2]) of the zero-range model must be recovered, we see that the lattice model's regular part tends to 
the zero-range model's regular part in the zero-range limit. 



9 



3. Finite-range continuous-space model 



Such models are used in numerical few-body correlated Gaussian and many-body fixed-node Monte-Carlo calcula- 
tions (see e. g. [El. l67l. 1801 - 183} and refs. therein). They are also relevant to neutron matter [84j . The Hamiltonian 
reads 

JV t N 

H = H + J2 E V ^ ( 21 ) 

i=l j=N f +l 

H being defined by ^ where A ri now stands for the usual Laplacian, and V(r) is an interaction potential between 
particles of opposite spin, which vanishes for r > b or at least decays quickly enough for r 3> b. The two-body 
zero-energy scattering state 4>(r) is again defined by the Schrodinger equation —(h 2 /m)A r cf) + V(r)cf> — and the 
boundary condition (|ll!12p . The zero-range regime is again reached for k typ b <C 1 with fc t yp the typical relative 
wavevector |l33l |. Equation (|20j) again holds in the zero-range regime, where A now simply stands for the zero-range 
model's regular part. 



B. Relations in the zero-range limit 



1. First order derivative of the energy with respect to the scattering length 



We now derive relations for the zero-range model. For some of the derivations we will use a lattice model and take 
the zero-range limit in the end. 
Three dimensions: 

Let us consider a wavefunction vb\ satisfying the contact condition ([]} for a scattering length a\. We denote by 
the regular part of vb\ appearing in the contact condition ((T|) . Similarly, yb 2 satisfies the contact condition for a 
scattering length a 2 and a regular part A? . Then, as shown in Appendix [B] the following lemma holds: 

ffV>2> - (ffVi,V<2) = — (- ~ -) (AW,AW) (22) 
where the scalar product between regular parts is defined by 

(AW,AW)=J2j ( 1] d d r k ) J ^.^(Ry,^)^)^? CH-tf.foW- (23) 

i<j k^i.j 

We then apply (|2"2"]l to the case where ipi and ip2 are TV-body eigenstates of energy E\ and E%. The left hand side of 
([22]) then reduces to {E 2 — Ei)(ipi\ip 2 ). Taking the limit a 2 — > a\ gives the final result 

dE iirti 2 . . .. 

-(A A) (24) 



d{— l/a) m 

for any eigenstate. This result is contained in the work of Tan [6ll l62j jl34| . Note that, here and in what follows, we 
have assumed that the wavefunction is normalized: (iftlip) = f . 
Two dimensions: 

The 2D version of the lemma (1221) is 



tyuHfo)- (Hik.,ih) = — ln(a 2 / ai ) (A^,A^), (25) 

m 



as shown in Appendix [Bj As in 31?, we deduce from the lemma the final result 

dE 2-wh 2 



, (1 , (A A). (26) 

d(ln a) m 
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2. Large-momentum tail of the momentum distribution 
The momentum distribution is denned in second quantization by 



rv(k) = (fit (k)c CT (k)} (27) 



where Co-(k) annihilates a particle of spin a in the plane-wave state |k) defined by (r|k) = e lk r . This corresponds to 
the normalization 



f d d k . Ar 



In first quantization, 



2 

d r 4 e- lk - r ^(ri,...,r A r)' 



(28) 



(29) 



where the sum is taken over all particles of spin o~, i.e. i runs from 1, to N-^ for a =t and from A^ + 1 to A" for a 
Three dimensions: 

The key point is that in the large-fc limit, the Fourier transform with respect to is dominated by the contribution 
of the short-distance divergence coming from the contact condition (TT]) : 

/ d? n e- ik ' ri V(ri, . / a?n e -* r - V —A^irj, (r k ) k ^). (30) 

J fe->oo J TYj 

From A(l/r) = — 4nS(r), we have the identity 
so that 

/ d 3 r ie - jk - r ^(r!,...,r w ) ~ ^ V e^A^r,, (r^,,). (32) 

Inserting this into (|29j) and expanding the modulus squared, the cross terms vanish in the large-fc limit, so that 

C = (4n) 2 (A,A) (33) 
where C = lim^oo fc 4 n CT (k). This can be rewritten using (|24p as: 



_ 4irm dE 

C V d(-l/o)' 1 j 



in agreement with Tan [6z 
Two dimensions: 



The 2D contact condition ([2]) now gives 

J d 2 r l e- ik - r ^(r 1 ,...,r N ) ^ J d 2 n ,e~ ik ^ £ ln(r y )A ij (r j , (r,)^). (35) 



From A(lnr) = 27ro"(r), we have the identity 



2tt 

d 2 re- ikr lnr = - TT , (36) 



so that 

| d 2 r, e-* r ^(r 1; . . . ,r N ) ^ £ e^^fa, (r,)^ij). (37) 

As in 3D this leads to 

C = (2n) 2 (A,A) (38) 

where C = limfc-^oo fc 4 n CT (k), and thus from ([24l) : 
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3. Short- distance asymptotic behavior of the pair distribution function 



The pair distribution function, Rivin g the probability density of finding a spin-f particle at point R + r/2 and a 
spin- 1 particle at point R — r/2, reads |135| 



. N r N 

g\f (R+£,R-£) = d d r 1 ...d d r N \ij(r 1 ,...,r N )\ 2 Y / E 6 (r + ~ r 4 ) 5 (r - | - r, 

i=l j=JV t +l 



JV T at 



e e /(n^)K ri '-> r * =R 



»=i i=iv t +i" fe^ij 
In what follows we consider the spatially integrated pair distribution function 



r r 

rj = K- r N 



G«(r) S /^4J(R+I, 



K -5 



(40) 



(41) 



(42) 



Three dimensions: 

Replacing the wavefunction in (|41j) by its asymptotic behavior given by the contact condition ([1]) immediately yields: 



Expressing (A, A) in terms of C trough (j3"3")l finally gives: 



r (2) fri r 



C 1 



(43) 



(44) 



In a measurement of all particle positions, the total number of pairs of particles of opposite spin which are separated 
by a distance smaller than s is 



so that from (|44 



^V P air(s)= / d*rG&(r) 



C 



^ pair (s) 



(45) 



(46) 



as obtained in [6JJ, [62| . 

Two dimensions: 

The contact condition @ similarly leads to 



G t?( r ) ~ 7^2 ln2r ' 

^ r^O (27r) 2 



After integration over the region r < s this gives 

ATpair(s) 



s->0 47r' 



2 i 2 

s In s. 



(47) 



(48) 



^. Expression of the energy in terms of the momentum distribution 



Three dimensions: 

As shown by Tan [61( , the total energy of any eigenstate has a simple expression in terms of the momentum distribution: 



rap Airma z -— ' 



ri 3 fc ft 2 fc 2 
(2tt) 3 2m 



C 

n CT (k) - 



(49) 
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or equivalently 



E - Strap = lim 

A— >oo 



H 2 C (\ 2A 



1 Ll\ \ ^— \ 
n 7r / -^-^ 



where 



and 



47TTO V a 7T 



C = lim fc Uo-(k), 

A:— >oo 



d 3 fc fi 2 fc 2 

fc<A (2tt)3 2m 



n CT (k) 



A 



E. 



trap 



E^ 1 



(50) 



(51) 



(52) 



is the trapping potential energy. A simple rederivation of this result is obtained using the lattice model (defined in 

Sec. Ill A"2|) : As shown in Section Hi C 31 one easily obtains Eq. (|100[) , which yields (|49p in the zero-range limit since 

D -> R 3 and e k -> h 2 k 2 /(2m) for b -> 0. 

Two dimensions: 

The 2L> version of (JSOJl is 



trap 



lim 

A— >oo 



n 2 C* , / aAe 7 
■ In 



27TTO 



E 



fc<A (2tt) 2 2m 



jv(k) 



(53) 



as was shown (for a homogeneous system) in 68]. This can easily be rewritten in the following forms, which resem- 
ble gH): 



E — Strap — — 



h 2 C , /age" 



In 



E 



d 2 fc fi 2 fc 



2i.2 r 



2irm *"\ 2 I ' ^ J (2tt) 2 2m 
where the Heaviside function 6 ensures that the integral converges at small k, or equivalently 



for any q > 0, 



E — Strap — — 



h 2 C , /age" 



hi 



2?rm V 2 



E 



2i.2 r 



d 2 k h 2 k 
(2tt) 2 2m 



i CT (k) 



C 



fc 2 (fc 2 +g 2 ) 



for any q > 0. 



(54) 



(55) 



To derive this we again use the lattice model. We note that, if the limit q — > is replaced by the limit b — > taken for 
fixed a, Eq. (Tl4l remains true (see App. |A"|); repeating the reasoning of Section Hi C 31 then shows that (I101[) remains 
true; taking the limit b — > finally gives 



E-R 



trap 



h 2 C , / age 1 
■In 



27rm V 2 

where g > is arbitrary; this can be rewritten as ([53 



E 



d 2 k h 2 k 2 

(2tt) 2 2m 



n„(k) - 



1 



k 2 k 2 — q 2 



(56) 



5. Short-distance expansion of the one-body density matrix 

The one-body density matrix is defined as 

^(R-|,R+0=(^(R-0^(k+0) (57) 

where Vv( r ) annihilates a particle of spin a at point r. Let us define a spatially integrated one-body density matrix 

G«(r) ee | d d i?. 9 W ( R - £,R + £) . (58) 

This is related to the momentum distribution by Fourier transformation: 

d d k 



G - (r) = /(2^ e "^ (k) - 



(59) 
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Three dimensions: 
As shown below, 

G«(r) = N„-£-r + 0(r 2 ), (60) 

and moreover the expansion can pushed to second order if one sums over spin and averages over three orthogonal 
directions of r: 

where the Ui's are three orthogonal unit vectors. This last relation, as well as its 2D version (|65|) . also hold if one 
averages over all directions of r uniformly on the unit sp here or unit circle. They generalize the result obtained in ID 
in 



551 ] . but the derivation is different from the ID case [136j . 



To derive (I60I61P we rewrite (|59"j) as 



OSS W + <«** - F + / <«"" ~ ') ("•« - F 

The first integral equals —{C/Sir)r. In the second integral, we use 



(62) 



e ,:kr -l = zkr- (k o r)2 +o(r 2 ). (63) 

The first term of this expansion gives a contribution to the integral proportional to the total momentum of the gas, 
which vanishes since the eigenfunctions are real. The second term is 0(r 2 ), which gives (ISO)) . Eq. (|6"Tj) follows from 
the fact that the contribution of the second term, after averaging over the directions of r, is given by the integral of 
k 2 [n a {k) — C/k 4 ], which is related to the total energy by (HH]). 
Two dimensions: 
As shown below, 

G«(r) = N<T + ^r 2 lnr + 0(r 2 ), (64) 
and for any pair of orthogonal unit vectors (ui,U2) 

2 



, , N+—r 2 

2—1 17 













+ 32 



. v . 2 (£-£ tl . ap )r 2 + (r 2 ) (65) 



(66) 



where 

T = 2 F 3 (l, 1; 2, 3, 3; -£) = -16 g (-!)' = 0.98625471 .... 

To derive (|64|65[) we rewrite (|59|) as 

G«(r) = ^ + /(r) + J(r) (67) 

with 

7 M=/|| (e ik --l)pfl(fc- ? ) (68) 

and 

J(r) = / (^ (glkr " 1} (" ff(k) " ~ 9) ) (69) 

where q > is arbitrary and the Heaviside function 8 ensures that the integrals converge. 
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To evaluate J(r) we use standard manipulations to get 



J(r) = ^r 2 [ln(qr) + 41] + 0(r 4 ) (70) 

87T 



1= / dx^t^, (71) 



where 



Jo being a Bessel function. Evaluating this integral with Maple gives 

7 - 1 - In 2 T 

X = - + (72) 

4 128 K ' 

where T is the hypergeometric function defined in ()66[) . 

Finally we evaluate J(r) using the same procedure as in 3D: expanding the exponential [see (|63p] yields an integral 
which can be related to the total energy thanks to (|54"1) . 

6. Second order derivative of the energy with respect to the scattering length 

We denote by \ip n ) an orthonormal basis of TV-body eigenstates which vary smoothly with 1/a, and by E n the 
corresponding eigenenergies. We will show that 

1 #E n (4nh 2 \ 2 _ p("'U(">)| 2 . 

w = l^rj JC^ m 3D (73) 



i ^ E l(^'U^))l 2 in2C (74) 



2d(lna) 2 \ m J ^ E n - E, 

n',E n ,^E n 

where the sum is taken on all values of n! such that E n i ^ E n . This implies that for the ground state energy Eq, 

3^ < ° " ^ (76) 

Eq.([75|) was intuitively expected 85]: Eq. (|46p shows that dEo/d(—l/a) is proportional to the probability of finding 
two particles very close to each other, and it is natural that this probability decreases when one goes from the BEC 
limit (—1/a — > — oo) to the BCS limit (—1/a — » +oo), i.e. when the interactions become less attractive |l37j . Eq. ([7rj)) 
also agrees with intuition [l38j j . 

For the derivation, it is convenient to use the lattice model (defined in Sec. Ill A"2|) : As shown in Sec III C 4l one easily 
obtains (|104l) . from which the result is deduced as follows. |0(O)| 2 is eliminated using (|17I18[) . Then, in 3D, one uses 

d 2 E n _ d 2 E n ( dg \ 2 , dE n d 2 g 



d(-l/a) 2 dg 2 \d(-l/a)J dg d{-l/a) 2 

where the second term equals 2go dE n /d{— 1/a) m/(4:irh 2 ) and thus vanishes in the zero-range limit. Similarly, in 2D 
one uses the fact that 



d 2 E n d 2 E n ( dg Q 



d(lna) 2 dg 2 \d(lna) 



(78) 



in the zero-range limit. 
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7. Time derivative of the energy 

We now consider the case where the scattering length ait) and the trapping potential U(r, t) are varied with time. 
The time-dependent version of the zero-range model (see e.g. (86[) is given by Schrodinger's equation 

d 

ih—ip(r 1) ...,r N ;t)=H(t)ii(r 1 ,...,r N ;t) (79) 
when all particle positions are distinct, with 



N r 



-A Tt +U(r t ,t) 
Zm 



(80) 



and by the contact condition (JlJ in 31? or in 2D for the scattering length a — a(t). One then has the relations 

fC " * 
dt Anm dt 



If = |^^^ + w*)iX>^*m*)> (si) 

»=1 

dE h 2 C dilii a) , , . . . „ . „ ^ . 

~rr = 7> ^K^I^ftPfri.tJWt)) in2D, (82) 



dt 2nm dt 



where E(t) — {^{t)\H{t)\vb{t)) is the total energy and E trsi p(t) = {ip(t)\ z~2iLi U(yi,t)\yb(t)) is the trapping potential 
energy. The relation (|5T|) was first obtained by Tan [62| . A very simple derivation of these relations using the lattice 
model is given in Section III C 51 Here we give a derivation within the zero-range model. 
Three dimensions: 

We first note that the generalization of the lemma (|2"2"|) to the case of two Hamiltonians Hi and if 2 with corresponding 
trapping potentials Ui(r) and ^(r) reads: 

(i>i,H 2 ^ 2 ) - (Hrfuth) = ^-(---) (AM,AW) + fa\ V [^(ri.t) - I7i(r is t)] |^ 2 ). (83) 

Applying this relation for = |t/>(£)) and ^2) = IVK* + St)) [and correspondingly ai = a(t), a 2 = a(t + ft) and 
Fx = #(*), H 2 =H(t + St)] gives: 

TO \ a W + "V / 

JV 

+<^(t)| 53 l/(ri,t + ft) -U-(ri,t)|^(t + **)). (84) 

i=l 

Dividing by St, taking the limit St — > 0, and using the expression (|33l) of [A, A) in terms of C , the right-hand-side of 
(|84")l reduces to the right-hand-side of (jSTj) . The left-hand-side of (|54"|) can be rewritten using Schrodinger's equation 
as iHd{ip(t)\tp(t + St)) /dt. Using Schrodinger's equation again to Taylor expand \if>(t + St)) in this last expression 
finally gives the result (|8Tj) . 
Two dimensions: 

The relation (|82p is derived similarly from the lemma 

2 fi 2 N 

(Vh,iW - Wi.^a) = — In(oa/oi)(^ (1) .^ (a) ) + V [C/ 2 (r 4 ,t) - l/k (*,*)] |^ 2 ). (85) 

m ' — ' 

i=l 



Relations for the lattice model 



In this Section, as well as in Section HTDl it will prove convenient to define C by 



„ 47TTO dE 

„ 27TTO dE 

c ^^^ inm (87) 
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This new definition of C coincides with the identities (|33|) and ((38]) of Section HlBl in the zero-range limit, as follows 
from (J2U2BJ). 

We will use the following lemma: For any wavefunctions ip, ip' , 

(i>'\W\iP) = \m\ 2 (A', A) (88) 

where A and A' are the regular parts related to ip and ip' through (119|) , and the scalar product between regular parts 
is naturally defined as the discrete version of 



The lemma simply follows from 

(^IW)=E E & (Ar - a) W)(ri,... J r < =r i> ...,r 3V ..,r JV ). (90) 

1, First order derivative of the energy with respect to the scattering length 
The Hellmann-Feynman theorem gives 

^ - = MW>. (91) 

ago ago 

But lemma (|88p with ip' = ip writes 

(ip\W\ip) = \4>{0)\ 2 (A, A). (92) 

Using the expressions (|1TI18|) of |0(O)| 2 , we conclude that 

dE Airh 2 . . . , . 

(A, A) in 3D (93) 

(A, A) in 2D. (94) 



d(—l/a) m 
dE 2-Kh 2 
d(lna) m 



2. Interaction energy 

The left-hand-side of (|92[) is obviously equal to the mean interaction energy E int divided by go; in the right-hand-side 
of ([55]). (A, A) can be expressed in terms of C using (|93I94I) and the definition (I86I87P of C: 

(A, A) = j£L i n 3D (95) 

(A A) = in 2D. (96) 

This gives 

n i \ ( f>(0)\ 2 ™3D (97) 

7 |0(O)| 2 in 2D. (98) 



Eint 


C 


go 


(4tt) : 


Eint 


C 


go 


(2tt) : 


:ives 





'h 2 \ l 

E iat = C I — — (99) 
>/ .go 



both in 3D and 2D. 
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3. Relation between energy, momentum distribution and C 



Here we show that 

E - Strap : 



E-E, 



trap 



h 2 c 



E 



d 3 k 

7 



in 3D 



,. H C . I aqe 
hm m 

«->0 lixm 



E 



d 2 fc 

(2^ 



n CT (k)-C- V 



2m€k 2m(ek — £q) 



(100) 

in 2D. (101) 



To derive this we start from the expression (|M|) of the interaction energy and eliminate 1/go thanks to (I13I14|) . The 
desired quantity E — -Etrap = E- lnt + Skin is then obtained from 



-Ekir 



E 



d d k 

{2ir) d 



e k n CT (k). 



(102) 



4- Second order derivative of the energy with respect to the coupling constant 

We denote by \ip n ) an orthonormal basis of iV-body eigenstates which vary smoothly with go, and by E„ the 
corresponding eigenenergies. We apply second order perturbation theory to determine how an eigenenergy varies for 
an infinitesimal change of go. This gives: 



E 



i<vviWn)r 



1 d 2 E n 

2 dg 2 ' ^ E n - E n , 
where the sum is taken over all values of n' such that E n i ^ E n . Lemma (l88l) then yields: 

2 



(103) 



1 d 2 E ^ /4( n ^P 

\i^ = ^ m E — 



(104) 



n',E n i^B„ 



5. Time derivative of the energy 
Equations (|81I82[) remain exact for the lattice model. Indeed, the Hellmann-Feynman theorem gives 



JY 



f-(?)-*m+(Ewh,». 



dt 



df 



(105) 



The result then follows by using the lemma (|88|) . the expressions (| 171181) of |</>(0)| 2 , and the expressions (|95I96[) of 
(A, A) in terms of C . 



6. On-site pair distribution function 
We define the spatially integrated pair distribution function 

G®W = 5.*- 5 



R 



(106) 



Using (|97I98I) and expressing the interaction energy in terms of gl, thanks to the second-quantized form (J9j) yields: 



G n(°) = 7£*rW°)l 2 in3I? 



(47T) 2 

C 



■|0(O)| 2 in 2D. 



(107) 
(108) 
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7. Pair distribution function at short distances 

The last result can be generalized to finite but small r, assuming that we are in the zero-range regime kt yp b <C 1 
(introduced at the end of Sec. Ill A 2|) : 



Gg)(r) 



C 



r«l/fc typ (27r) ; 



■\4>{y)\ 2 in 2D. 



(109) 
(110) 



Indeed, the expression (j4Tj) of gj^ in terms of the wavefunction is valid for the lattice model with the obvious 
replacement of the integrals by sums, so that 



N t N 



R i=l j=/V t +l (r fc )**i,i 



) • • • )M 



r r 

r, = R +-,..., Yj = R - 2 ' ■ • ■ ' rN 



(111) 



For 7- <C 1/fctyp, we can replace ^ by the short-distance expression (|20l) . assuming that the multiple sum is dominated 
by the configurations where all the distances (r^ — R| and ruu 1 are much larger than b: 



G$(t)~(A,A) Wr)! 1 



(112) 



Expressing (A, A) in terms of C thanks to (|95I96|) gives (jl09lll0|) . 



8. Momentum distribution at large momenta 
Assuming again that we are in the zero-range regime k iyv b <C 1, we will show that 



Mk) ~ C 



(—)' 

V 2me k / 



both in 3D and in 2D. We start from 



v(k)=E £ 



(JV-l) 



E^ kr ^(ri,...,r w ) 



(113) 



(114) 



We are interested in the limit k 3> &typ- Since ^(ri, . . . , rjv) is a function of rj which varies on the scale of l/fc typ , 
except when is close to another particle Yj where it varies on the scale of b, we can replace ip by its short-distance 
form PTJ): 



E^e- <k -*V(ri, ..., rjv ) ~ £(k) E < ' k , '-M»V faW, 



(115) 



where 0(k) = (r 1 0) = 6 d e~ ,k ' r </>(r). Here we excluded the configurations where more than two particles are at 
distances < b, which are expected to have a negligible contribution to (|114l) . Inserting (|115p into (|114l) . expanding 
the modulus squared, and neglecting the cross-product terms in the limit k S> fctypj we obtain 



u CT (k)~|0(k)| 2 (4,A). 



(116) 



Finally, </>(k) is easily computed for the lattice model: for k ^ 0, the two-body Schrodinger equation (|A1|) directly 
gives 0(k) = -5o</>(0)/(2e k ), and 0(0) is given by (|15I16I) . which yields (fTT5j) . 



D. Relations for a finite-range interaction in continuous space 



We recall that in this Section, C is again defined by (|86I87|) . 
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1. Interaction energy 
As for the lattice model, we find that the interaction energy is proportional to C: 



E\ nl = / ,/ : YlV)Mr) 1 in :\D (]1.7) 



£fat = J d 2 rV{r)\ct>{r)\ 2 in 2D. (118) 

It was shown in [65j that this relation is asymptotically valid in the zero-range limit in 31?. Here we show that it 
remains exact for any finite value of the range and we generalize it to 2D. 
For the derivation, we set 

V(r) = g W(r) (119) 

where go is a dimensionless coupling constant which allows to tune a. The Hellmann-Feynman theorem then gives 
Eiat — godE/dgo. The result then follows by writing dE/dgo = dE/d(—l/a) ■ d(—l/a)/dgo in 31? and dE/dgo = 
dE/d(lna) ■ d(ln a)/dgo in 2D, and by using the definition (|86I87[) of C as well as the following lemmas: 

s °^r = I d2rV{rmr)f m2D - (121) 

To derive these lemmas, we consider two values of the scattering length ai, i — 1,2, and the corresponding 
scattering states <pi and coupling constants go,i- The corresponding two-particle relative- motion Hamiltonians are 
Hi = —(h 2 /m) A r + g 0yi W(r). Since H^i = 0, we have 



R— >oo 



lim / d?r (faHifc ~ faHifa) = 0. (122) 



r<R 



The cont ribution of the kinetic energies can be computed from Ostrogradsky's theorem and the large-distance form 
of <j> |l39| . The contribution of the potential energies is proportional to go. 2 — ffo,i- Taking the limit 02 — > a\ gives 
the results (|120I121I) . Lemma (|120l) was also used in [65| and the above derivation is essentially identical to the one 
of [65| . For this 31? lemma, there also exists an alternative derivation based on the two-body problem in a large 
box [140| . 

2. Relation between energy and momentum distribution 

Three dimensions: 



h 2 C ^ f d 3 k h 2 k 2 

E — Strap — "J h / 



4irma ^ J (27r) 3 2m 
where 4>'(k) = 0(fc) + a - 1 (27r) 3 (5(k) = J d 3 r e^'VM with 



n - (k) -(^^' (k)|2 



(123) 



0'(r) = <j>(r) + -. (124) 



a 

This is simply obtained by adding the kinetic energy to (|117j) and by using the lemma: 

o , , vl0 47rft 2 f d 3 k h 2 k 2 -,, , 

d 3 rV r )\4> r 2 = / ? — 3 \4>' [k 2 . 125 

mo J (27r) m 

To derive this lemma, we start from Schrodinger's equation — (ft 2 /m)A0 + T^(r)</> = 0, which implies 

^2 



J d 3 rV{r)\cj){r)\ 2 = ^ J <* 3 



rcj)A(j). (126) 
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On the other hand, applying Ostrogradsky's theorem over the sphere of radius R, using the asymptotic expression 
(fTTj) of <j) and taking the limit R — > oo yields 

ri 3 r0A0=-^- J d 3 r(V</>) 2 . (127) 

We then replace by V</>'. Applying the Parseval-Plancherel relation to di(f>, and using the fact that <j>'(r) vanishes 
at infinity, we get: 

d?r{V^>) 2 = / -^k 2 \4>>(k)\ 2 ( 128 ) 

The desired result (I125[) follows. 

Two dimensions: 



h 2 C , /R\ f d 2 k h 2 k 2 



n.(k) - ^|^(fc)| 2 



(129) 



_r.-s.oo I 27rm la/ ^— ' / (2tt) 2 2m 
where 4>' R (k) = J d 3 r e~* kr <f>' R {r) with 

W = W) - HR/a)] 9(R - r). (130) 
This follows from (]118|) and from the lemma: 

d 2 r^ r Mr 2 = lim In - - / —— <t>' R (k 2 . 131 

fl^oo [ m \ a J J \2ttY m J 

The derivation of this lemma again starts with (|126[) . Ostrogradsky's theorem then gives [139] 



-R— s-oo V a 



d 2 r4>Ac))= lim {27rln - - / tf r (V<£) ^ ■ (132) 



r<i? 



We can then replace J r<R d 2 r (V0) 2 by J" <i 2 r(V<// fl ) 2 , since <fi' R {r) is continuous at r = -R [139] so that V0^, does 
not contain any delta distribution. The Parseval-Plancherel relation can be applied to di<f>' R , since this function is 
square- integrable. Then, using the fact that 4>' R (r) vanishes at infinity, we get 



and the lemma (] 1 3 1 [) follows. 



<l'rlY»> ; r I ^k 2 \4>' R (k)\ 2 , (133) 



3. Pair distribution function at short distances 



In the zero-range limit kt yp b <C 1, the short-distance behavior of the pair distribution function is given by the same 
expressions (|109I110|) as for the lattice model. Indeed, Eq. (|112p is derived in the same way as for the lattice model; 
one can then use the zero-range model's expressions (|33|38|) of (A, A) in terms of C, since the finite range model's 
quantities C and A tend to the zero-range model's ones in the zero-range limit. 

4- Momentum distribution at large momenta 
In the zero-range regime kt yp b -C 1 the momentum distribution at large momenta k 3> fctyp is given by 



n„ 



77, 



(k) ~ J^\m\ 2 in 3t9 (134) 
C 



(k) * J^y\^)\ ^ 2t9. (135) 
Indeed, Eq. (|116[) is derived as for the lattice model, and (A, A) can be expressed in terms of C as in Subsec lTlD 31 
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E. Derivative of the energy with respect to the effective range 

We now show that, in 3D, the leading order finite-range correction to the zero-range model's spectrum is given by 
the model-independent expression 



I7) = 27r E f d3R f{H d 3 r k )A^{v k ) k ^ 3 ) 



N 



1 = 1 



E-— A,!-!-— > A,., - > ru-,i A zj (R,(r k ) k ^) 

(136) 

where r e is the effective range of the interaction potential, the derivative is taken in r e — 0, the function A is assumed 
to be real without loss of generality, and in the sum over I we have set = Tj — R. To obtain this result we use a 
modified version of the zero-range model, where the boundary condition ([1} is replaced by 



V>(ri,...,rjv) = 



I 771 

- + 7^2? r e ) Aij (Rij, {r k )k^,j) + 0(nj), 



(137) 



where 



£ = E-2U(H ij ) - [ u ( r ") 



Aij (R*ij, (rfe)fc^i, - 



Am 



2m 



A,, :R,j.{r k ) kj ,,,}. (138) 



Equations (I137I138[) generalize the ones already used for 3 bosons in free space in [83, 1881 (the predictions of [87J and 
[88j have been confirmed using different approaches, see [89| and Refs. therein, and [90, [Sll respectively). Such a 
model was also used in the two-body case, see e.g. [921-94] . and the modified scalar product that makes it hermitian 
was constructed in [95]. 

For the derivation of (I136|) . we consider an eigenstate ip\ of the zero-range model, satisfying the boundary condition 
(fTj) with a scattering length a and a regular part A^\ and the corresponding finite-range eigenstate ip2 satisfying 
(|137|138j) with the same scattering length a and a regular part A^ 2 K As in App. |5]we get (|B3j) . as well as (|B6j) with 
1/ai - l/a 2 replaced by m£r e /(2h 2 ). This yields (jTB^l) . 



F. Generalization to statistical mixtures and to thermodynamic equilibrium in the canonical ensemble 

The above results, summarized in Tables [HI Hill and IIV1 hold for any eigenstate (apart from lines 8-9 of Tab. [TT] 
and line 11 of Tab. IIII[) . Thus they can be generalized straightforwardly to statistical mixtures of eigenstates. The 
relation for the time derivative of E (Tab. [IT] line 9) holds for any time-evolving pure state, and thus also for any 
time-evolving statistical mixture. 

We turn to the case of thermal equilibrium in the canonical ensemble. We shall use the notation 

A^{7 1/a in3D (139) 
1 ilna in 2D. 



a. First order derivative of E. The thermal average in the canonical ensemble dE/dX can be rewritten in the 
following more familiar way, as detailed in Appendix [C] 



~d\) ~ \ dXJ T ~ \ dx) t 



(140) 



where (...) is the canonical thermal average, F is the free energy, U — E is the mean energy and S is the entropy. 
Taking the thermal average of (|24I26[) thus gives 

dF ) J^U^M (141) 



dX J j, \dX J 

The other results (I44I47I34I39I49I53I) are generalized to finite temperatures in the same way. 



22 



b. Second order derivative of E . Taking a thermal average of the expression (173I74P we get after a simple ma- 
nipulation: 



l(d 2 E\ (ATth 2 y 1 



2 V d\ 2 J \ m 
where Z — ^ n exp(—f3E n ). This implies 



2Z 



E 



n,n';E n jtE n 



E n — E n i 



.|(A( n '),AW)| s 



d 2 E\ 

n» <0 - 



Moreover one can check that 



which implies 



f cPF 



(d 2 E 
\~d~X 2 



f d 2 F 
\~dX 2 



dE 
~dX 



< 0. 



dE 
~dX 



<0, 



(142) 



(143) 



(144) 



(145) 



In usual cold atom experiments, however, there is no thermal reservoir imposing a fixed temperature to th e gas, one 
rather can achieve adiabatic transformations by a slow variation of the scattering length of the gas [96Ml00j . One also 
more directly accesses the mean energy U of the gas rather than its free energy, even if the entropy is also measurable 

35L The second order derivative of U with respect to A for a fixed entropy is thus the relevant quantity to consider 

101]. As shown in the appendix [Cl one has in the canonical ensemble: 



d 2 U 
d\ 2 



_ fd 2 E\ 



[Gov(E, 4g )] 2 - Var(£)Var( 



dE \ 
d\ ) 



k B TV&r(E) 



(146) 



where Var(X) and Cov(X, Y) stand for the variance of the quantity X and the covariance of the quantities X and 
Y in the canonical ensemble, respectively. From the Cauchy-Schwarz inequality [Cov(X, Y)] 2 < Var(X)Var(F), and 
from the inequality (|143l) . we thus conclude that 



d 2 U\ n 



(147) 



To be complete, we also consider the process where A is varied so slowly th at th ere is adiabaticity in the many-body 
quantum mechanical sense: The adiabatic theorem of quantum mechanics [102] implies that in the limit where A is 
changed infinitely slowly, the occupation proba bilitie s of each eigenspace of the many-body Hamiltonian do not change 
with time, even in presence of level crossings [l03j |. We note that this may require macroscopically long evolution 
times for a large system. For an initial equilibrium state in the canonical ensemble, the mean energy then varies with 
A as 



-p a E n {\ ) 



z 



E n (X) 



(148) 



where the subscript refers to the initial state. Taking the second order derivative of (I148[) with respect to A in 
A = Aq gives 



72 pquant 
" -^adiab 

d\ 2 



( d2E \ n 



(149) 



Finally, we compare the result of isentropic transformation (| 1461) to the one of the adiabatic transformation in the 
quantum sense (|149p . They differ by the second term in the right hand side of (|146[) . A priori this term is extensive, 
and thus not negligible compared to the first term. We have explicitly checked this expectation for the Bogoliubov 
model Hamiltonian of a weakly interacting Bose gas. The discrepancy between (|146p and p49p indicates that the 
limit of infinitely slow transformation does not commute with the thermodynamic limit. More explicitly, we see that 
if A is varied so slowly that the quantum adiabaticity (|148p is achieved, one cannot assume any more that the system 
follows a sequence of thermal equilibrium states with a constant entropy; in practice this may require evolution times 
which grow exponentially with the system size. 
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III. SPINLESS BOSONS 



The wavefunction ip( r i> ■ • • i t n) is now completely symmetric. 

Our results for bosons are shown in Table fVl An obvious difference with the fermionic case is that there are no more 
spin indices in t he p air distribution function g^> and in the momentum distribution n(k). Accordingly, Eqs. (I40I28[) 



are replaced by [1411 



and 



7 «(r+E R-E) = J dv 1 ...dv N \i){v 1 ,...,v N )\ 2 ^5(K+ r --v^5(n- T --v J ) (150) 



"'%n(k) = N. (151) 



(27T) 



An important difference with the fermionic case is that in 3D, the Efimov effect occurs 71|, and the zero-range model 
is defined not only by the contact condition ([T]) and the Schrodinger equation but also by a boundary condition 
in the limit where three particles approach each other: There exists a function B, hereafter called three-body regular 
part, such that 



^ (ri "--' rjv) ^oi sin 



l s ol ln ^7 



*(n)B(C,r 4> ...,rjv). (152) 



where Rt is the so-called three-body parameter and is an additional parameter of the zero-range model; C = (ri + 
r 2 + r 3)/3 is the center of ma ss of particles 1,2 and 3; R and f2 are the hyperradius and the hyperangles associated 
with particles 1,2 and 3 [l42j | We recall the definition of R and f2: The Jacobi coordinates are defined by r = T2 — ri 



and p — (2r 3 — ri — r 2 )/-\/3; then, R = \J (r 2 + P 2 )/2, and tt = (a,r,p) with a = arctan(r//o), r = r/r and 
p = p/p. so — i ■ 1.00624. . . is Efimov's transcendental number, it is the imaginary solution of (|192[) . $(ri) is the 
normalized hyperangular part of an Efimov states' wavefunction: <P(£l) — 4> So (O) / y/ (0s o |0 So ) [142] . We recall that, in 
the present case (bosons with zero total angular momentum), cf> So (&l) = (1 + Q) sin [so (f — a)] /[V^n sin(2a)] where 
Q = P13 + P23 and Pij exchanges particles i and j. The hyperangular scalar product is defined by (4>\<j>) = J dfl |^(r2)| 2 

with J dtt = 2 J^ 2 da sin 2 (2a) J df J dp, where dr and dp are the differential solid angles. Its value is given in App. [D] 
For N = 3 it is well established that this model is self-adjoint and that it is the zero-range limit of finite-range 
models, see e.g. [75| and references therein. A recent numerical study of several finite range models predict ed fo r 
N ~ 4 the existence of tetramers of e nergi es weakly depending on the model for a fixed th ree-body parameter [104| , 
an existence confirmed experimentally jl05| | . This numerical study, together with other ones [l06l I lOTl j claim that there 
is no need to introduce a four-body parameter in the zero-range limit, implying that the here considered zero-range 
model is self-adjoint for N — 4. Here we consider an arbitrary value of N such that the model is self-adjoint. 

The main relations for zero-range interactions are displayed in Table |Vj Moreover we note that the relations for 
finite-range interactions, given in Tables Hill and HVl for fcrmions, can be easily generalized to the bosonic case. 



A. Relations which are analogous to the fermionic case 

The derivations of all relations of Table [V] are completely analogous to the fermionic case, except for lines 4 and 6 
of the left column (3D case). 

The first result in Table [V] was first obtained in [75| in the case N = 3. A simple way to derive it for any N is to 
use the lattice model and to apply the reasoning of Scc lII C ll The key point is that in the limit of a lattice spacing b 
much s malle r than |o|, the three-body parameter corresponding to the lattice model is equal to a numerical constant 
times b [143| . Thus, varying the coupling constant go while keeping b fixed is equivalent to varying a while keeping 
R t fixed: 

dE _( dE \ d(-l/a) 



d 9o \d(-l/a) J Rt dg 

where the lattice mod el's (dE/d{— l/a))n t tends to the zero-range model's one if one takes the zero-range limit while 
keeping R t fixed [144| . The same reasoning explains why the second derivative in the last line of Table [V] also has 
to be taken for a fixed Rt- 
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B. Derivative of the energy with respect to the three-body parameter 

The Efimov effect also gives rise to the following new relation between the derivative of the energy with respect to 
the three-body parameter and the three-body regular part: 



dE \ h 2 V3 



dlnR, 



= ^^\s \ 2 N(N-l)(N-2) JdcJ dv i ...dr N \B(C,T4,...,r N )\' 



(154) 



This is simil ar to the relation (IM1) between the derivative with respect to the scattering length and the (two-body) 
regular part [l45| . We will first derive this relation using the zero-range model in the case N = 3, and then using a 
lattice model for any N. 

1. Derivation using the zero-range model for three particles 

We consider two wavefunctions ipi, ip2, satisfying the two-body boundary condition ([T]) with the same scattering 
length a, and satisfying the three-body boundary condition with different three-body parameters Rn, Rti- The 
corresponding three-body regular parts are denoted by B±, B 2 . We show in the App. [El that 



(ih,Hih) - (HM 2 ) = ^^r^sin 

m lb 



Nm — 

tUx 



j ' dCBl(C)B 2 (C), (155) 



which yields (|154j) by choosing ipi as an eigenstate of energy Ei and taking the limit Rti — > Rti- We note that ipi and 
ip2 do not satisfy lemma l|22p because they are too singular for R — > 0. 

2. Derivation using a lattice model 

We now rederive (|154|) using a lattice model which is analogous to the model defined in Sec III A 2l except that the 
Hamiltonian now contains a three-body interaction term: 

/d 3 k 
j (2^ £kct ( k ) c W + E b 3 U(r)(yj^)(r) + g £ 6 3 (^W^)(r) + h £ & 3 (^V^,t^)( r ). (156) 

The scattering length is related to go as in Sec lII A2l We define the zero-energy three-body scattering state 
(/>o(ri, r 2 , r 3 ) as the solution of i?|</>o) = for a = oo, with the boundary condition 



/f>o(ri,r 2 ,r 3 ) ~ sin 



|so|ln — 
tit 



(157) 



in the limit where all interparticle distances tend to infinity. This defines the three-body parameter R t (b,ho) for 
the lattice model. In order to derive (11541) for any desired values of a and R t , we first choose b in such a way that 
Rt{b, ho = 0) is equal to the desired R t divided by e n7r /' s °', the zero-range limit n — > oo being taken in the end. We 
then choose go to reproduce the desired value of a. We then change ho to a small non-zero value, keeping fixed b and 
go (and thus also a). The Hellman-Feynman theorem writes: 

^ = ^i 3 ((^W#^(r))=iV(]V-l)(JV-2) b^ N -^\iP(r,r,r,r 4 ,...,r N )\ 2 . (158) 

r r 4 ,...,r\v 

For the lattice model we define the three-body regular part B through: 

ip(r, r, r, r 4 , . . . , r N ) = <fo(0, 0, 0) B(r, r 4 , . . . , r N ); (159) 

in the zero-range limit, we expect that this lattice model's regular part tends to the regular part of the zero-range 
model defined in (|152[) , as was discussed for the two-body regular part A in Sec. Ill A 21 We thus have, in the zero-range 
limit: 

dE \ 2 f dh 



i)(]uIh yj =N(N-l)(N-2)\M0,0,0)\ 2 ( g^y ) / d 3 rd 3 n...d 3 r N \B(r,r 4 ,...,r N )\ 2 . (160) 

It remains to evaluate the expression between brackets: This is achieved by applying (|160[) to the case of an Efimov 
trimer in free space, where the regular part can be deduced from the expression for the normalized wavefunction given 
in App. |Dj 
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C. Momentum distribution for an Efimov trimer 



For an Efimov trimer state at rest, for an infinite scattering length, we show in Appendix [F] that the atomic 
momentum distribution has the asymptotic expansion 



n(k) 



C_ 



D 

fc5 



2|s |ln(fcV3/ 



where sq 



1.00624... solves 



s cos 



(aw/2) - 8/%/3sin(s7r/6) = 0, 



(161) 



(162) 



the energy of the trimer E tl i m = —h 2 K^/m depends on the three-body parameter R t as specified in (|D1[) . and the 
quantities C, D and ip are derived in the appendix [F] The crucial point is that D ^ 0: The momentum distribution 
has a slowly decaying oscillatory subleading tail. 

The calculations performed in appendix [F] also allow a straightforward numerical calculation of the atomic mo- 
mentum distribution for an Efimov trimer, both for low values of k, see Figfjji, and for high values of fc, see FigOJ) 
showing how n(k) approaches the asymptotic behavior (|16ip . We have also derived in appendix [F] the exact value in 
k = 0: 



n(k = 0) 



55.43379775608.. 



(163) 




100 




FIG. 1: For a free space Efimov trimer at rest composed of three bosonic atoms of mass m interacting via a zero range, 
infinite scattering length potential, atomic momentum distribution n(k) as a function of k. (a) Numerical calculation from the 
expressions derived in the appendix[F] (b) Numerical calculation (solid line) and asymptotic behavior (| 161 p (dashed line), with 
the horizontal axis in log scale. The unit of momentum is Ko, such that the trimer energy is —tf/rriKQ. 



D. Breakdown of the energy-momentum relation in the zero-range model 

1. A non-converging integral 

As a consequence of (I16ip . the integral / d 3 k k 2 [n(k) — C/k 4 } is not well-defined: After the change of variables 
x = In k, the integrand behaves for x — > oo as a linear superp osition of e l ^ x and e ~ l|s ° |:E , that is as a periodic 
function of x oscillating around zero. This was overlooked in [681 ]. 



2. Failure of a naive regularisation 



At first sight, however, this does not look too serious: one often argues, when one faces the integral of such an 
oscillating function of zero mean, that the oscillations at infinity simply average to zero. More precisely, let us define 
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the cut-off dependent energy of the Efimov trimer (here 1/a = 0) 

d 3 k H 2 k 2 
(2tt) 3 2m 



E(A) = f 

Jk<A 



(164) 



For A — > oo, E(A) is asymptotically an oscillating function of the logarithm of A, oscillating around a mean value E. 
The naive expectation would be that the trimer energy -Btrim equals E. This naive expectation is equivalent to the 
usual trick used to regularize oscillating integrals, consisting here in introducing a convergence factor e _1,ln (' c / Ko ) in 
the integral without momentum cut-off and then taking the limit r\ — > + : 



lim 

?7->0+ 



d 3 k h 2 k 2 
(2tt) 3 2m 



n(k) 



-rj ln(fc//<o) 



E. 



(165) 



We have decided to really test this naive regularisation. As we now show, remarkably, it actually fails to give the 
correct energy of the Efimov trimer. 

a. Numerical evidence for E ^= E. We first performed a numerical calculation of E(A), using the results of 
Appendix IF1 to perform a very accurate numerical calculation of n(k). The result is shown as a solid line in Fig|2] 
We also developed a more direct technique allowing a numerical calculation of E(A) without the knowledge of n(fc), 
see Appendix [G] The corresponding results are represented as + symbols in Fig|5] and are in perfect agreement with 
the solid line. As expected, E(A) is asymptotically an oscillating function of the logarithm of A, oscillating around a 
mean value E. This may be formalized as follows. We introduce an arbitrary, non-zero value fc m in of the momentum, 
and we define 



Sn(k) = n(k) 
6n(k) = n(k) 



^ for k <C k m m 
C_ D_ 

k I + ¥ 



2 1 s 1 ln(fcV / 3/Ko) + V 



for k > k„ 



(166) 
(167) 



The introduction of k mm ensures that the integral of k 2 5n(k) over k converges around k = 0. The subtraction of the 
asymptotic behavior of n(k) up to order 0(1/ k 5 ) for k > k m i n ensures that the integral of k 2 Sn(k) over R 3 converges 
at infinity. As a consequence we get for A > k min the splitting 
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With the change of variable x = ln(fc\/3/Ko)j one can calculate the second integral explicitly. Since the first integral 
converges in the limit A — > +00 we obtain 



E(A) =E- 
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From this last equation and the numerical calculations of n(k) first up to k = l OOO/ tp and then up to A; 
get two slightly different values of E which give an estimate with an error bar 146| : 



(169) 

(170) 
5500/to> we 



E ~ 0.89397(3)£ t rim- (171) 

The key point is that E significantly differs from E^i^: the naive regularisation does not give the correct value of the 
trimer energy! 

b. Physical explanation and expression of E. We now show that an analytical expression for E may be obtained, 
which confirms the numerical conclusion and has the crucial advantage of explaining the physics behind the discrepancy 

E 7^ -Etrim- 

The first step is to realize what happens in a finite range interaction model, when one takes the zero range limit. E.g. 
in the lattice model, for a non-zero value of the lattice spacing b. one readily realizes that an exact energy-momentum 
relation holds for an arbitrary number of bosons even in three dimensions. Repeating the reasoning of Sec. Ill B 41 
while keeping finite the lattice spacing b, we get 
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(172) 
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ln(A/K ) 

FIG. 2: Cut-off dependent energy 15(A) as denned in (|164|l for a free space infinite scattering length Efimov trimer with a 
zero range interaction, as a function of the logarithm of the momentum cut-off A. Solid line: numerical result obtained via a 
calculation of the momentum distribution n(k). Symbols +: direct numerical calculation of 15(A) as exposed in the Appendix 
[Gl Dashed sinusoidal line: asymptotic oscillatory behavior of 15(A) for large A, obtained in omitting 0(1/A) in (|169[) . Dashed 
horizontal line: mean value E around which 15(A) oscillates at large A. The values of E obtained analytically (|181[) and 
numerically (|171|) are indistinguishable at the scale of the figure, and clearly deviate from the dotted line giving the true energy 
Etiim of the trimer, exemplifying the failure of a at a first sight convincing application of an energy-momentum relation for 
bosons in three dimensions. The unit of momentum kq is such that the true trimer energy is 15 tr i m = —T^K^/m. 



where 

and the prime denotes quantities calculated within the lattice model. If one then takes the zero-range limit by 
repeatedly dividing b by the discrete scaling factor [144], so as to ensure a well-defined value for the three-body 
parameter R t , the lattice model's quantities E' , E' t , C and n'(k) converge to the zero-range model's quantities E, 
Etrap, C and n(k); however one cannot simply remove the primes and replace the integration domain by R 3 in (|172l) . 
because this would lead to an ill-defined integral, as we have seen. This paradox is due to the fact that the finite 
range interaction b in the lattice model has two effects, that conspire to ensure that the energy-momentum relation 
(|172l) is correct: (i) it introduces a momentum cut-off of order 1/6, and (ii) it changes the large momentum tails of 
the momentum distribution in a way that tends to zero when b — > for fixed k, and that still can have a non-zero 
impact on the energy-momentum relation since the integral of k 2 /k 5 is UV divergent in 3D. The previous failure of 
the energy-momentum relation for the zero range model is thus due to the fact that we have introduced a cut-off (|164|l 
or a regularisation (|165l) in the integral over k without subsequently modifying in an appropriate way the coefficient 
of the 1/fc 5 part of n{k). 

In a second step, we introduce a consistent way of performing a UV regularisation, that will both show up in the 
integral over k in the energy-momentum relation and affect the tail of the momentum distribution. To be explicit, we 
turn back to the example of N — 3 bosons forming an Efimov trimer in free space, of energy i/trim = ~fr 2 KQ/m, for an 
infinite scattering length. In the zero range model, when the positions i"i and r 2 of the particles 1 and 2 symmetrically 
converge to a common point R 12 in space, the three-body wavefunction diverges as U(2|r 3 — R 12 |/v / 3)/(— 47rr 12 ) where 
the function B is known, see Appendix [F] The one body momentum distribution has an explicit integral expression 
in terms of the Fourier transform B(k) of this function B, as shown in that Appendix. The idea is then to introduce 
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a smooth regularisation simply replacing B(k) with [1471 ] 

B v (k) = B(k) e - I ' ln [^ /1+fe2/K » +fe/K °] (174) 

and eventually take the limit 77 —> + . Performing the replacement (|174l) in the expression (|F15IF16IF17IF18IF19[) of 
the momentum distribution leads to a modified momentum distribution n n (k), with an asymptotic behavior n^ symp (re) 
for k — s- +00 that is modified as compared to (|161j) : n v (k) = n^ symp (k) + 0(1/ k 6 ) with 

k— >+oo 

< ymp (fc) = + p |Aj + A, cos [2|s | ln(fcV3/ Ko ) + p„J } • (175) 

In the limit 77 — + one has to recover (|16ip so that C v —> C, D v —> D, ip v —> ip and D v — > 0. The expressions of C v 
and D v are given in the Appendix U and confirm these requirements. What shall play a crucial role in what follows 
is that, however, D v is not zero for 77 > 0, it vanishes linearly with rj in 77 = 0. For this consistent regularisation, the 
energy-momentum relation holds in the limit of vanishing 77 for the Efimov trimer, as we prove in Appendix IHl 
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The definitions (|166p . (|167l) are then modified as 



n v( k ) - -gk 



(176) 



C_ 

k 4 

5n n (k) = n v (k) - n^ ymp (k) for k > fc min . (178) 



6n n (k) = n v (k) - tj for fc < /c min (177) 



This results in the splitting 



f ri 3 k h 2 k 2 h 2 r + °° 

E v= 7o~T3 ~o — Sn ri (k) + - — 2 — / dxe- 2 ^[D n + D. n cos(2\ So \x + ^)] (179) 

where the change of variable x — ln(fc-\/3/ kq) was used so that £ m in = \n(k m i n V3 / kq) ■ For 77 — > + , we can replace in 
the right hand side of (I179P Sn v (k) with 5n(k) since the first integral converges absolutely, but we cannot exchange 
the limit and the integration in the second integral. After explicit calculation of this second integral, we take r\ —> + 
and we recognize E from (|170|) so that 

E trim = E+-^- lim ^L. (180) 

As detailed in the Appendix HI D v and the above limit may be expressed as single integrals, which allows to evaluate 
E with a good precision: 

E = 0.8939667780883 . . .-E trim (181) 
This confirms the numerical estimate (|171l) . 



IV. ARBITRARY MIXTURE 



In this Section we consider a mixture of bosonic and/or fermionic atoms with an arbitrary number of spin com- 
ponents. The N particles are thus divided into groups, each group corresponding to a given chemical species and to 
a given spin state. We label these groups by an integer a € {1, . . . Assuming that there are no spin-changing 
collisions, the number N a of atoms in each group is fixed, and one can consider that particle i belongs to the group a 
if i € I a , where the I^s are a fixed partition of {1, ... , N} which can be chosen arbitrarily. For example, a possible 
choice is I\ — {1, . . . ,Ni}; I2 = {N\ + l,...,Ni + N2}; etc. The wavefunction ?/>(ri, . . . ,rpj) is then symmetric 
(resp. antisymmetric) with respect to the exchange of two particles belonging to the same group I a of bosonic (resp. 
fermionic) particles. Each atom has a mass mi and is subject to a trapping potential £/i(rj), and the scattering length 
between atoms i and j is ay . We set = m a and ay = a aa ' for i 6 l a and j G . The reduced masses are 
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(i aa > — m a m a i j{m a + m a t). We shall denote by P aa i the set of all pairs of particles with one particle in group a and 
the other one in group a', each pair being counted only once: 



P aa , = € {h X la') U (I a . Xl a ) / l< j} . 



(182) 



The definition of the zero-range model is modified as follows: In the contact conditions (|l|2j) the scattering length 
a is replaced by , and the limit r%j — > is taken for a fixed center of mass position R.y = (m^r^ + mj-rj)/(mi + m.,-); 
moreover Schrodinger's equation becomes 



JV r 



E 



2m. 



•Ar, +Ui(Ti) 



4' = Eip. 



Our results are summarized in Table IVI[ where we introduced the notation 



(i,j)eP„ 



[] d d r k ) J d d %^g ) (R ij -,(r fe )^ J OMg ) (R ii ,(r fe )^ iJ ). 



(183) 



(184) 



We note that since = a ff ' ff there are only n{n + l)/2 independent scattering length, and the partial derivatives 
with respect to one of these independent scattering lengths are taken while keeping fixed the other independent 
scattering lengths. 

In 31? the Efimov effect can occur, e.g. if the mixture contains a bosonic group, or at least three fermionic groups, or 
two fermionic groups with a mass ratio strictly larger than a critical value 13.6 . . . [H3 |. In this case, as in the previous 
Section, the derivatives with respect to any scattering length have to be taken for fixed three-body parameter (s), and 
the relation between E and the momentum distribution (line 4 of Table IVI[) breaks down [l48| . This relation was 
first obtained in [68[ in 3D, and in 2D for Fermi-Fermi mixtures. Here we express the first sum in a more explicit 
way in terms of the partial derivative of the energy, and we point out the breakdown of this relation in presence of 
the Efimov effect. 

The derivations are analogous to the ones of Sections |TT] and IIIII The lemmas (I22I25P are replaced by 



2irh 2 { 1 



irh 2 



l 

a 2 



(4 (1 U (2) )^ in 3D 



ln(a 2 / ai )(A«,A( 2 >) c 



(185) 



in 2D. 



The pair distribution function is now defined by 

fl£ 2 J/(u,v) =/ dr x ...dr N |V>(ri,..., rjv)| 2 ^ S (u - r^) 6 (v - 

The Hamiltonian of the lattice model used in some of the derivations now reads 

W, 



a<a' 



(186) 



(187) 



where 



N 



i=l 



2m 



■A ri + Ui(ri) 



with (r|A r |k) = -A: 2 (r|k) and 



b~ d . 



(188) 



(189) 



In the formulas of Sec. Inland App. |A"1 coming from the two-body scattering problem, one has to replace go by go,a-a-', 
a by a aa i and m by 2[i r7r7 i. Denoting the corresponding scattering state by 4> aa i(v), the lemma (|88[) becomes 



(^n>\W aa ,\4 n ) = |<^w(0)| 2 (A^'\A^) C 



(190) 
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V. APPLICATIONS 

In this Section, we apply some of the above relations, first to the three-body problem and then to the many-body 
problem. We consider the unitary limit a — oo in three dimensions. 

A. Three-body problem: corrections to exactly solvable cases and comparison with numerics 

Here we use the known analytical expressions for the three-body wavefunctions to compute the corrections to the 
spectrum to first order in the inverse scattering length 1/a and in the effective range r e . 

1. Universal eigenstates in a 



The prob lem of three identical spinless bosons [1081 . Il09j or two-component fermions (say = 2 and N± = 
1) [lOSL \llo\ is exactly solvable in the unitary limit in an isotropic harmonic trap U(r) = 1/2 moj 2 r 2 . Here we restrict 
to zero total angular momentum, with a center of mass in its ground state, so that the normali zatio n constants of the 
wavefunctions are also known analytically (75|. Moreover we restrict to universal eigenstates |l49{ . The spectrum is 
then given by 

E = E cm + (s + 1 + 2q)Hui (191) 
where E cm is the energy of the center of mass motion, s belongs to the infinite set of real positive solutions of 

_SC0S ( S f) +7] ~^ S ' m ( S ^) = ° (192) 

with r\ = +2 for bosons and —1 for fermions, and q is a non- negative integer quantum number describing the degree 
of excitation of an exactly decoupled bosonic breathing mode [86llllll ]. We restrict for simplicity to states with q = 0. 

a. Derivative of the energy with respect to 1/a. Injecting the expression of the regular part A of the normalized 
wavefunction [75[ into the relation (|24l) or its bosonic version ( Table fVl line 1) we obtain 



,)E r( S + l/2)V2 S sin( S |) [ml 



8{-l/a) 



T(s + 1) - cos (sf) + sf sin (sf ) + rj^ cos (sf 



For the lowest fermionic state, this gives (8E / 8(1 / a)) a =oo — — 1.1980 y^h 3 u>/m, in agreement with the value —1.19(2) 
which we extracted from the numerical solution of a finite-range model presented in Fig. 4a of [8lj |. where the error 
bar comes from our simple way of extracting the derivative from the numerical data of [8l| . 

b. Derivative of the energy with respect to the effective range. Using relation (|136|) . which holds not only for 
fermions but also for bosonic universal states, we obtain 

(^-) = = r (B -l/2)^-l/2)sin(W2) JJ—^ m 

\8r e J a r(s + l)2\/2 [-cos(stt/2) + svr/2 • sin(svr/2) + n2Tr/(3V3) ■ cos(s?r/6)] 

For bosons, this result was derived previously using the method o f 18711 and found to agree with the numerical solution 
of a finite-range separable potential model for the lowest state [7a|. For fermions, (|194[) agr ees w ith the numerical 
data from Fig. 3 of [81] to ~ 0.3% for the two lowest states and 5% for the third lowest state [150]; (|194|) also agrees 
to 3% with the numerical data from p. 21 of [75| for the lowest state of a finite-range separable potential model. All 
these deviations are compatible with the estimated numerical accuracy. 

2. Derivative of the energy of an Efimov trimer with respect to 1/a. 

The same relation ( Table IVl line 1) can be applied to Efimov trimers in free space. Using the expression of the 
normalized three-body wavefunction (see App. [D| we get 

8E \ = (_tf E \ 1/2 ^sinh(| S o|7r/2)tanh(| So |7r) 



0(-1AO/h, V "» J CO sh(| S o|7r/2) + -^sinh(| So |7r/2)-^cosh(| So |7r/6) 



2 011111 VI°oi»/^ -75' 

h 2 \ 1/2 

= 2.11267159347... x [ E) . (196) 



FIG. 3: Pair distribution function g!^(r) = 
(V'|(r)V)|(0)V'4(0)i/>^(r)) of the homogeneous unitary 
gas at zero temper atur e. Crosses: fixed-node Monte-Carlo 
results from Ref. |ll3| . Line: analytic expression (|200|l . 
where the value £ = 0.95 was taken to fit the Monte-Carlo 
results. The arrow indicates the range b of the square-well 
interaction potential. 



FIG. 4: One-body density matrix g£)(r) = {4>t (r)V> CT (0)) 
of the homogeneous unitary gas at zero temperature: com- 
paris on b etween the fixed-node Monte-Carlo results from 
Ref. 114] (continuous curve) and the analytic expres- 
sion (|199|) for the small-fcj?r expansion of g^J up to first 
order (dashed straight line, red online) and second order 
(dotted parabola, blue online) where we took the value 
£ = 0.95 extracted from the Monte-Carlo data for gW , see 
Fig. El 



This confirms and refines the value of the coefficient 2.11 estimated numerically in [73[. We note that this numerical 
coefficient is simply A'(— 7r/2)/| sp|, w here A(£) is Efimov's universal function that was estimated numerically in [73[ 
and computed very precisely in |112j . Our analytical calculation gives the exact value of the derivative 

A'(-tt/2) = 2.125850069373..., (197) 

to be compared with the numerical estimate A'(— n/2) ~ 2.12 in (73[. 

The expression (|195l) can also be obtained from the tail of the momentum distribution: The expression of the 
coefficient C follows from Eq. flFll) . and (|195[) is recovered using the relation on the second line of the left column of 
Table El 



B. Unitary Fermi gas: comparison with fixed-node Monte-Carlo 



For the homogeneous unitary gas (i.e. the two-component Fermi gas in 31? with a = oo) at zero temperature, we 
can compare our analytical expressions for the short-distance behavior of the one-body density matrix gjj} and the 
pair distribution function g!i) to the fixed- node Monte-Carlo results published by the Trento group in [82l . Ill3l . Ill4| . 

In this case, g£r(R — r/2, R + r/2) and <?1?(R — r/2, R + r/2) depend only on r and not on a, R and the direction 
of r. Expanding the energy to first order in l/(kpa) around the unitary limit yields: 



E = E £ 



kpa 



(198) 



where Eq is the ground state energy of the ideal gas, £ and £ are universal dimensionless numbers, and the Fermi 
wavevector is related to the density through kp = {i^n) 1 / 3 . Expressing C in terms of ( thanks to (|34I198[) and 
inserting this into (161[) gives 



l-^k F r--^(k F r) 2 + 
10 10 v ; 



For a finite interaction range b, this expression is valid for b <C r -C kp 1 15l|. Equation (|109[) yields 



(199) 



(200) 
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The interaction potential used in the Monte-Carlo simulations [§2, Ill3l . 1 1 1 4j ] is a square- well 



v{r) = " (I ) ^ ^ r < b (2Q1) 



for r > b 



for which the zero-energy scattering state is 



4 > {r) = \ sin (^)A forr<6 
I 1 /r for r > b 

and the range 6 was taken such that nb 3 = 1CT 6 i.e. k F b = 0.0309367 .... Thus we can assume that we are in the 
zero-range limit k F b <C 1, so that (|199I200|) are applicable. 

Figure [3] shows that the expression (|200|) for g& fits well the Monte-Carlo data of [l!3| | if one adjusts the value of 
C to 0.95. This value is close to the value £ ~ 1.0 extracted from (|198p and the -E(l/a)-data of [82]. 

Using £ = 0.95 we can compare the expression (|199l) for g^J with Monte-Carlo data of [114] without adjustable 
parameters. Figure 2] shows that the first order derivatives agree, while the second order derivatives are compatible 
within the statistical noise. This provides an interesting check of the numerical results, even though any wavefunction 

satisfying the contact condition ([1]) would lead to g£~} and g^ functions satisfying (|44I60I) with values of C compatible 
with each other. 



C. Dependence of the unitary gas critical temperature on the interaction range 

The key property underlying f| 136[) is that the leading order change of each eigenstate energy of a spin 1/2 Fermi 
gas due to a small but non-zero interaction range is linear in the effective range r e of the interaction potential, 
which a coefficient which is model independent. As a consequence, the leading order change of the thermodynamical 
potentials, and even of the critical temperature T c of the Fermi gas, are also expected to be linear in r e , with model 
independent coefficients. 

This expectation can be tested with the Quantum Monte Carlo data of [4l[ and [3, where the critical temperature 
of the unitary gas was calculated for two different, finite range interaction models [l52j . In Figj5l we plot the Monte 
Carlo data as a function of kpr e , where r e is the effective range of the corresponding model. Our expectation is that 
all the data lie on the same straight line for low enough A;i?|r e |, which is indeed essentially the case (considering the 
size of the error bars). We then conclude that the shift of T c due to the interaction range is of order 

ST 

—Z ~ 0.12fc F r e , (203) 

-LF 

at low r e and in a model independent way. This results in a relative shift at the percent level for typical experiments 
on lithium, where kpr e s» 0.01. 



VI. CONCLUSION 



We derived relations between various observables for N particles of arbitrary masses and statistics in an external 
potential with zero-range or short-range interactions, in continuous space or on a lattice, in two or three dimensions. 

Some of our results generalize the ones of [55l l6ll l62l l65l . f6Sl j : The large- momentum behavior of the momentum 
distribution, the short-distance behavior of the pair correlation function and of the one-body density matrix, the 
derivative of the energy with respect to the scattering length or to time, the norm of the regular part of the wavefunction 
(defined through the behavior of the wavefunction when two particles approach each other), and, in the case of finite- 
range interactions, the interaction energy, are all related to the same quantity C; and the difference between the total 
energy and the trapping potential energy is related to C and to a functional of the momentum distribution (which 
is also equal to the second order term in the short-distance expansion of the one-body density matrix) . For Efimov 
states with zero-range interactions, we found that this last relation breaks down, because the large-momentum tail of 
the momentum distribution contains a subleading oscillatory term. 

We also obtained entirely new relations: The second order derivative of the energy with respect to the inverse 
scattering length (or to the logarithm of the scattering length in two dimensions) is related to the regular part of 
the wavefunctions, and is negative at fixed entropy; the derivative of the energy of a universal state with respect to 



33 



0.25 



0.2 



0.15 



0.1 



0.05- 







1 1 1 


1 1 

s- 

y 




y 

S — 

* 


y 

'. I--''''' 


- 


•at 

Jiff 

1 


i 



-0.5 







0.5 



k„ r 

F e 



FIG. 5: Critical temperature T c of the unitary Fermi gas as a function of the effective range r e of the interaction potential, as 
given by the Quantum Monte Carlo results of [4l| for the Hubbard model (symbols with r e < 0) and of [H| for a continuous 
space model (symbols with r e > 0). The dashed line corresponds to a linear fit of the data over the interval kpr e £ [—0.8, 0.45]. 
Here UbTf = h 2 kp/(2m) is the Fermi energy of the ideal gas with the same density as the unitary gas. 



the effective range of the interaction potential is also related to the regular part; and the derivative of the energy 
of an efimovian state with respect to the three-body parameter is related to the a three-body analog of the regular 
part. Applications were presented in three dimensions for an infinite scattering length: the derivative of the energy 
with respect to the inverse scattering length was computed analytically and found to agree with numerics for Efimov 
trimers; the same was done for universal three-body states in a harmonic trap, not only for the derivative of the 
energy with respect to the inverse scattering length but also with respect to the effective range; existing fixed-node 
Monte-Carlo data for the unitary Fermi gas were checked to satisfy exact relations. Also the derivative of the critical 
temperature of the unitary gas with respect to the effective range, expected from our results to be model- independent, 
is estimated from the Quantum Monte Carlo results of |4lj ]. 

The relations obtained here may be used in various other contexts. For example, the result (11471) on t he si gn of the 
second order derivative of E at constant entropy is relevant to adiabatic ramp expe riments [35l I36l.l98l l99l . Ill5| , and the 
relation (I107P allows to directly compute C u sing diagrammatic Monte-Carlo Ill6i . C is directly related to the closed- 
channel fraction in a two-channel model ffriL l66j . which allows to extract it j66| from experi mental ph otoassociation 
measurements |32|. C also pl ays an important role in the theory of radiofrequency spectra (65l . Ill7l - ll21 1 and in finite-a 
virial theorems 163, Ii23.li23j. 

We can think of several generalizations of the relations presented here. All relations can be rederived in the case 
of periodic boundary conditions. The relations for finite-range models, obtained here for two-component fermions, 
can be generalized to arbitrary mixtures. The techniques used here can be applied to the one-dimensional case to 
generalize the relations of [55j . For two-channel or multi-channel models one may derive relations other than the ones 
of [64H66{ . In presence of the Efimov effect, the derivative of the energy with respect to the three-body parameter can 
easily be related to the short-distance behavior of the third order density correlation function thanks to (|154I152I) ; 
moreover the asymptotic behavior (|161[) of the momentum distribution is expected to hold for any state satisfying 
the three-body boundary condition (|152[) . with a coefficient D of the subleading tail in (jl61[) related to dE/d(\nR t ) 
through a simple proportionality factor. Indeed, the singularities of the wavefunction at short interparticle distances 
are generally related to large momentum tails. 



34 



Acknowledgments 

We thank S. Tan, M. Olshanii, E. Burovski, N. Prokof'ev and R. I gnat for usefu l disc ussions, as well as S. Giorgini, 
J. von Stecher and M. Th0gersen for sending numerical data from [8ll Ill3l . Ill4 Il24l |. F.W. is supported by NSF 
under Grant No. PHY-0653183. Y.C. is a member of IFRAF. 

Note 

While completing this work, we b ecam e awa re of unpublished notes by Tan where Eqs. (|39|82|) were obtained 
independently using the formalism of [l25| [l53j | . 

Appendix A: Two-body scattering for the lattice model 

For the lattice model defined in Sec. Ill A 2^ we recall that <fi(r) denotes the zero-energy two-body scattering state 
with the normalization (|llll2p . In this Appendix we derive the relation ()13|14[) between the coupling constant go 
and th e sca ttering length, as well as the expressions (|15I16I17I18I) of </>(0). Some of the calculation resemble the ones 
in [781.IT261]. 

We consider a low-energy scattering state $ q (r) of wavevector q <C b^ 1 and energy E = 2e q ~ H 2 q 2 /m, i.e. the 
solution of the two-body Schrodinger equation (with the center of mass at rest): 

(H + V)\^) =E\^) (Al) 

where Ho = J D d d kj(2ii) d 2ek|k)(k| and V — go\r = 0)(r = 0|, with the asymptotic behavior 

p iqr 

* q (r) = e lqr + / q + ... in 3D (A2) 



$ q (r) = e^ r -/ q W- e** + ... in 2D. (A3) 

r->oo \/ 8irqr 



Here / q is the scattering amplitude (in 2D the present definition is from |l27| . and y/i = e l7T ^ 4 ), which in the present 
case is independent of the direction of r as we will see. We then have the well-known expression 



where G — (E + i0 + — -ff) -1 , which, since 

is equivalent to 

where 



|$ q ) = (l + GV)|q) (A4) 
G = Go + GqVG, (A5) 
|$ q ) = (l + G T)|q) (A6) 



T = V + VGV. (A7) 

Indeed, (|A4I) clearly solves (|A1[) . and one can check [using the fact that (r|Go|r = 0) behaves for r — > oo as 

-m/{^h 2 ) e lqr /r in 3D and ~m/h 2 v / i/(8Trqr)e lqr in 2D] that satisfies (IA2IA3|) with 

/q = ~4S^ 3(r = °™ [U3D (A8) 

/ q = ^6 2 (r = 0|T|q) in 2D. (A9) 

Using (fATl) and (KE\i one gets 

1 f d d k 1 



<r = 0|T|q) = b 



-d 



go J D (27r) d E + i0+ - 2e k 



(A10) 
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In 3D the scattering length in defined by / q — > —a, which gives the relation (fl3|) between a and go. In 2D, 

/q _ ln( (/ ae^/2)-i7r/2 + o(l) (All) 



where a is by definition the 2D scattering length. Identifying the inverse of the right-hand-sides of Eqs. (|A9[) and 
(|Alip and taking the real part gives the desired (|14p . We note that Eqs. (|A11I14I) remain true if q — > is replaced by 
the limit b — > taken for fixed a. 
To derive (|15I16[) we start from 

V|$) =T(£ + iO + )|q) (A12) 
which directly follows from (|A4j) . Applying (r = 0| on the left and using (|A8IA9|) yields 

ffo<f q (0) = ~^-/q in 3D (A13) 
. 9o $ q (0) - ^-/ q in 2D. (A14) 
In 3D, we simply have <ft = — a~ 1 lim$ q , and the result (|T5|) follows. In 2D, the situation is a bit more tricky because 

q— >0 

lim$ q (0) = 0. We thus start with q > 0, and we will take the limit q — > later on. At finite q, we define </> q (r) as 

being proportional to $ q (r), and normalized by imposing the same condition (|12[) than at zero energy, but only for 
l)«r«f Inserting (|A11|) into (|A14I) gives an expression for $ q (0). To deduce the value of <j>(0), it remains to 
calculate the r- independent ratio <^ q (r)/$ q (r). But for r>fcwe can replace </> q (r) and < I > q (r) by their values within 
the zero-range model (since we also have b -C q^ 1 ) which we denote by (/> q R (r) and $ q R (r). The two-body Schrddinger 
equation 

fi 2 

A$™ = E Vr > (A15) 

implies that 

$ ZR = gl q.r + ^(1) ( qJ .j ( A1 g) 

where Af is a constant and i?o is an outgoing Hankel function. The contact condition 

3A/ $™(r) = AMr/a) + 0(r) (A17) 

together with the known short-r expansion of the Hankel function [l28| then gives 

A = — ^ -. (A18) 

ln(gaeT/2)-i7r/2 V ; 

Of course we also have <l> q R /</> q R = A, which gives (fl6|) . 

Finally, Eqs. (|17ll8p are obtained from (|15I16[) using the relations d(m / (An fi 2 a)) /d(l/g ) = 1 in 3D and 
d(l / go) / d(ln a) — —m/(2nh 2 ) in 2D, which are direct consequences of the relations (|13ll4p between go and a. 



Appendix B: Derivation of a lemma 

In this Appendix, we derive the lemma (|22p in three dimensions, as well as its two-dimensional version (|25p . 
Three dimensions: 
By definition we have 

h 2 r' N 

{^i,H^)-{H^x,^) = -— ^...dSv^I^Ar^-^Ar^il- (Bl) 
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Here the notation f ' means that the integral is restricted to the set where none of the particle positions coincide [l54| . 
We rewrite this as: 

j_2 N pi p 

(ih.,Hi> i )-(Hi> 1 ,^) = -—Y l (T\d 3 r k )lim d 3 n [^A^ - ^A r ^J] . (B2) 

lm i= i J V ' e_> ° J{nM¥*,r»>e} 

We note that this step is not trivial to justify mathematically. The order of integration has been changed and the 
limit e — > has been exchanged with the integral over r^. We expect that this is valid in the presently considered 
case of equal mass fermions, and more generally provided the wavefunctions are sufficiently regular in the limit where 
several particles tend to each other. 

Since the integrand is the divergence of ^*V ri ^2 — V>2V r4 ^>*, Ostrogradsky's theorem gives 

j-2 w pi 

{ll>l,Hip 2 )-{mi,^) = —Y, (U^jl^Yl § [^iVr^-lfeV, 4 ^]-dS (B3) 

where the surface integral is for belonging to the sphere S e (rj) of center Tj and radius e, and the vector area dS 
points out of the sphere. We then expand the integrand by using the contact condition, in the limit = e — > taken 
for fixed rj and fixed (Tk)k^i,j- Using Ry = Yj + eu/2 with u = (r^ — v^)jr^ we get 

£-£) A V + r^^' + Oic ) (B4) 

V„*„ = + i [V„„A<;" - u (u ■ V„„4")] + 0(1) (B5) 

where n equals 1 or 2, and the functions and Vr^.A^™'' are taken at (r^, (rk)k^n,j)- This simply gives 

§ mVr^2 ~ ^ 2 V ri Vi*] • dS = q 4tt (- - 4 1} *4 2) + 0(e) (B6) 

because the leading order term cancels and most angular integrals vanish. Inserting this into (|B3I) gives the desired 
lemma (1221 . 
Two dimensions: 



The derivation is analogous to the 3D case. In (|B3j) . the double integral on the sphere of course has to be replaced 
by a simple integral on the circle. Instead of (|B4IB5|) . we now obtain, from the 2D contact condition ©, 

i> n = ln(e/o») A n) + O(elne) (B7) 

e— »-0 J 



V ri V„ = -A^+OQne), (B8) 



which gives 

/ Wl^v^i ~ V^V, 4 Vi] • dS = 27rln(a 2 /a 1 )4l ) *^ 2) +0(eln 2 e) (B9) 
and yields the lemma (|25l) . 

Appendix C: First and second order isentropic derivatives of the mean energy in the canonical ensemble 

One considers a system with a Hamiltonian H(X) depending on some parameter A, and at thermal equilibrium 
in the canonical ensemble at temperature T, with a density operator p — exp(—(3H)/Z. In terms of the partition 
function Z(T, A) = Tie- f3H( - x \ with /3 = l/(fc B T), one has the usual relations for the free energy F, the mean energy 
U = Tv(pH) and the entropy S = — fc B Tr(plnp): 

F(T, A) = -k B TlnZ(T,\) (CI) 
F(T,A) = U(T,X)-TS{T,X) (C2) 
d T F(T,X) = -S(T,X). (C3) 
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One now varies A for a fixed entropy S. The temperature is thus a function T(A) of A such that 

5(T(A),A) =ct. (C4) 
The derivatives of the mean energy for fixed entropy are then: 

§), s Iptom)] (o») 
tsf), s ^m^i. <C6) 

Writing (|C2I) for T = T(A) and taking the first order and the second order derivatives of the resulting equation with 
respect to A, one finds 



(dU 



= d\F(T(\), \) (C7) 

s 



*^ =^ ( r W ,A)-^m!. (C8) 



dAVs 9|F(T(A),A) 

It remains to use (|C1[) to obtain a microscopic expression of the above partial derivatives of F, from the partition 
function expressed as a sum Z = ^ n e~P En over the eigenenergies n of the Hamiltonian: 

dxF(T,X) = ^ (C9) 



^,r,A) . g-^v„(f ) (cio) 

* F(r - A » = (C11 » 

. ,ci2) 



The expectation value (. . .) stands for a sum over the eigenenergies with the canonical probability weights, and Var 
and Cov are the corresponding variance and covariance. E.g. E = U and 



d 2 E d 2 E n e-^ 

n 

Cov(E,dE/dX) = Y,E n ^^—-E^. (C14) 



dA Z dA 

n 

Insertion of JCSJ into ([C3) gives ([140]) . Insertion of (|C10IC1 LIC12I) into (Cg) gives (Tl46| 



Appendix D: Normalized wavefunction of an Efimov trimer 

In this Appendix we recall the wavefunction of an Efimov trimer and give the expression of its normalization 
constant. We consider an Efimov trimer state for three spinless bosons of mass m interacting via a zero range infinite 
scattering length potential. In order to avoid formal normalisability problems, we imagine that the Efimov trimer is 
trapped in an arbitrarily weak harmonic potential, that is with a ground state harmonic oscillator length eto arbitrarily 

h 2 K 2 

larger than the trimer size. In this cas e, th e energy of the trimer is essentially the free space energy -Etrim = — ^f-, 
k > 0. According to Efimov's theory |l29l | 

Ko = ^ e - ? /| S ol e Argr(i+ So )/M (Dl) 
Rt 

where R t > is a length known as the three-body parameter, the quantum number q may take all values in Z and 
the purely imaginary number so = i|so| is such that 

|s | cosh(| So k/2) = A sinh(|s k/6), (D2) 
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so that | sq | = 1.00623782510 . . . The corresponding three-body wavefunction ^ may be written as 



*(ri,r 2 ,r 3 ) Vcm(C) [^(n 2 , |2r 3 - (n + r 2 )|/V3) + ^(r 23 , |2n - (r 2 + r 3 )|/VS) + ^(r 3 i, |2r 2 - (r 3 + ri)|/V3)] , 

(D3) 

where C = (ri + r 2 + r 3 )/3 is the center of mass position of the three particles and the parameterization of ip is 
related to the Jacobi coordinates r = r 2 ri and p = [2r 3 — (ri + r 2 )]/v / 3- In our expression of vp, -0cm is the 
Gaussian wavefunction of the single particle ground state in the harmonic trap, normalized to uni ty, and ip is a 
Faddeev component of the free space trimer wavefunction. The explicit expression of ip is known [l29j : 



ip(r, p) = 



K S0 (K y/r 2 + p 2 ) sin [s (| - a)] 



4tt 



(r 2 +p 2 )/2 



sin(2a) 



(D4) 



where K SQ is a Bessel function and a = atan(r/p). The normalization factor ensuring that ||'J / || 2 = 1 may 
be calculated explicitly: One first performs the change of variables (ri,r 2 ,r 3 ) —> (C,r,p), whose Jacobian is 
D(ri,r 2 ,r 3 )/D(C, p,r) = (— 1/3/2) 3 . To integrate over r and p one introduc es hy perspherical coordinates in which 
the wavefunction separates; one then faces known integrals on the hyperradius [l30j | and on the hyperangles [83] • This 
leads to (zl: 




3tt 2 



2kq cosh(|so|7r/2) 



cosh(|so|7r/2) 



sinh(|soK/2) 



47T 

3Vl 



cosh(|s k/6) 



We also recalled, as promised in the main text, the value of the hyperangular scalar product derived in [7£ 



sol0s o ) = — sin(s 7r/2) 
so 



7r Air 

cos(so7r/2) — so— sin(so7r/2) -= cos(so7r/6) 

2 3v 3 



(D5) 



(D6) 



Appendix E: A lemma for three bosons in the zero-range model 

Here we prove the relation (|155[) . The first step is to express the Hamiltonian in hyperspherical coordinates [142] 

h 2 fVs' 3 



(Vi,i^2)-(i^i> 2 ) = -A-f^j J^dRR 5 JdnJdC 



tf(^ + 52s + § + H*-rc^} (E1) 

— |y dRR 5 J dnAc(R,Cl) + J dfldCA R {n,C) 

dRR 5 J dCA u (R,C)^ (E2) 



where 



A C (R, n) = J dC U* iA c ^ 2 - Wi ^ H 
A R (n,C) = J dRR 5 U*i ( " 



(E3) 

V>2 - W{ H } (E4) 



\dR 2 RdR 

An(R,C) = Jdfl j^^Va- V-2^1 , (E5) 

To being a differential operator acting on the hyperangles and called Laplacian on the hypersphere. 

Clearly Ac(R, fi) = 5 JdC Vc • {ipiVcfa — ^Vc^i} = 0, since the ipiS are regular functions of C for every 
(R, ft) except on a set of measure zero. 

In what follows we will use the following simple lemma: if <&i(R) and <I > 2 (i?) are functions which decay quickly at 
infinity and have no singularity except maybe at R = 0, then 

f d ««° {*; [m + Ira) - « M - - fc R {^H - ^} (E6 » 
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where ^(R) = R 2 
We now show that 



An(R, C) = for any C and R > 0. (E7) 



We will use the fact that ip\ and 1^2 satisfy the two-body boundary condition with the same a, and apply lemma (|22l) . 
More precisely, we will show that for any smooth function f{R, C) which vanishes in a neighborhood of R = 0, 

J dRR 5 J dCf(R,C) 2 An(R,C) = 0; (E8) 

this clearly implies (|E7[) . To show (|E8j) we note that 




y di?i? 5 y dC/(i?,C) 2 A^,C) = -|J- J dRR 5 J dn J dc|(M)§(/^a)-[^^^]}, 

(E9) 



which can be rewritten as 



1 d 3 r 1 d 3 r 2 d 3 r 3 {(M*)i?(/V'2) - [Vi ^ *h]} +^ J dRB? J dtt J dC 



{jn?-Tim + l Ac] (EII,! 



The first integral in this expression vanishes, as a consequence of lemma (J22J). This lemma is indeed applicable to the 
wavefunctions ftpi: They vanish in a neighborhood of R = (see the discussion below (|B2[1 ). moreover they satisfy 
the two-body boundary condition for the same value of the scattering length a (as follows from the fact that R varies 
quadratically with r for small r). The second integral in (|E10[) vanishes as well: The contribution of the partial 
derivatives with respect to R vanishes as a consequence of lemma (|E6|) . and the contribution of Ac vanishes because 
the ftpi's are regular functions of C. 

Finally, Ar can be computed using lemma (|E6|) and the boundary condition (I152p . yielding (|155[) . 

Appendix F: First two terms of the large-fc expansion of the momentum distribution of an Efimov trimer 

Here we show that the momentum distribution of an Efimov bosonic trimer state of energy —h 2 K^/m (at rest and 
for an infinite scattering length) has the asymptotic expansion (|161[) with 

8-7T 2 sinh(|s |7r/2) tanh(|s |7r) 



C/n = , , ' U| ' ' v^LJ. = 53.09722846003081 . . . (Fl) 

cosh(|sok/2) + i^ s inh(| So |7r/2) - ^ cosh(| So k/6) 

D/kI ~ -89.26260 (F2) 

if ~ -0.8727976 (F3) 



1. Three-body state in momentum space 

We start from the three-body wavefunction in position space ^ given in Appendix [Dj To obtain the momentum 
distribution of the Efimov trimer, we need to evaluate the Fourier transformation of ^ . Rather than directly using 
(|D4p . we take advantage of the fact that the Faddeev component ip obeys Schrodinger's equation with a source 
term. With the change to Jacobi coordinates, the Laplacian operator in the coordinate space of dimension nine reads 
Ei=i A r, = | A c + 2 [A r + A p ] so that 

- [k 2 -Ar- A p ] i>(r, p) = S(r)B( P ). (F4) 
The source term in the right hand side originates from the fact that 

<Kr,p) ~ (F5) 
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for a fixed p, this 1 jr divergence coming from the replacement of the interaction potential by the Bethe-Peierls contact 
condition. Taking the Fourier transform of (|F4[) over r and p leads to 



<Kk,K) 



B(K) 



K + Kq 



(F6) 



where the Fourier transform is defined as B(K) = / d 3 pe- t ' K p B{p). B{p) is readily obtained from flD4]) by taking 
the limit r — » 0: 



- -^(47r) 1 / 2 isinh(| So |7r/2)^^. 
The Fourier transform of this expression is known, see relation 6.671(5) in [l30j, so that 



2tt 5 / 2 



^(if 2 + K 2 ) 1 ^ 



(K 2 + k 2 ) 1 / 2 + K 



(K 2 + /tg)V2 + A' 



K. 



(F7) 



(F8) 



What we shall need is the large K behavior of B(K). Expanding (|F8|) in powers of Ko/K gives 

B(K) = N^^- [(2K/ Ko r so ~ c.c] + 0(1/X 4 ). 



K 

When necessary one may further use the relation 



( Ko /2) so = (R t V2) 



-■o r(i + a ) 
|r(i + So ) 



(F9) 



(F10) 



that can be deduced from (|D1|) . 

The last step is to take the Fourier transform of (|D3[) . using the appropriate Jacobi coordinates for each Faddeev 
component (or simply by Fourier transforming the first Faddeev component using the coordinates (C,r, p) given 
above and by performing circular permutations on the particle labels). This gives 



*(ki,k 2 ,k3) = (^J ^ CM (ki+k 2 + k 3 ) [vl(|k 2 - kil/2, V3|k 3 -(k 1 +k 2 )/2|/3) 



^(|k 3 - k 2 |/2, V3|k! - (k 2 + k 3 )/2|/3) + ^(|ki - k 3 |/2, V3|k 2 - (k 3 + k 1 )/2|/3)l . (Fll) 



2. Formal expression of the momentum distribution 



To obtain the momentum distribution, it remains to integrate over k 3 and k 2 the modulus square of (|F11[) . In the 
limit Ko a ho +oo, one can set 



|^c M (ki + k 2 + k 3 )| 2 = (2^) 3 <5(k! + k 2 + k 3 ). 
Integration over k 3 is then straightforward: 



(F12) 



i(ki) = 3(V3/2) 



(2tt) 3 



^(|k 2 - ki|/2, V3|ki + k 2 |/2) + ^(|k 2 + ki/2|, y/Zkx/2) + ^(|k a + k 2 /2|, V3fc 2 /2) 



(P13) 

The factor 3 in the right hand side results from the fact that, in this article, we normalize the momentum distribution 
ri (k) to the total number of particles (rather than to unity). One further realizes that the sum of the squared moduli 
of the arguments of ip is constant and equal to k\ + k 2 + ki • k 2 for each term in the right hand side. One uses (|F6I) . 
thus putting the denominator in (|F6[) as a common denominator, to obtain 



l( k 1 ) = 3(V3/2)«/|| 



B(V3|ki + k 2 |/2) + S(V3*i/2) + B(V3k 2 /2) 



{k 2 + k 2 + ^ • k 2 + k 2 ) 



(F14) 
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For simplicity, we have assumed that the normalization factor is purely imaginary, so that B(K) is a real quantity. 

In the above writing of n(ki), the only "nasty" contribution is £?(\/3|ki + lt2 1/2) ; the other contributions are "nice" 
since they only depend on the moduli k\ and k 2 . Expanding the square in the numerator of (|F14|) . one gets six terms, 
three squared terms and three crossed terms. The change of variable k 2 = — (k 2 + ki) allows, in one of the squared 
term and in one of the crossed term, to transform a nasty term into a nice term. What remains is a nasty crossed 
term that cannot be turned into a nice one; in that term, as a compromise, one performs the change of variable 
k 2 = — (k 2 + ki/2). We finally obtain the momentum distribution as the sum of four contributions, 



n(k x ) = m(ki) + n//(ki) + nj W (k x ) + n/v(M (F15) 



with 



d 3 k 2 


B 2 {V3k 1 /2) 


(F16) 


(2tt) 3 


(fc 2 + fcf + ki • k 2 + k 2 ,) 2 


d 3 k 2 


2S 2 (73fc 2 /2) 


(F17) 


(2tt) 3 


(fc 2 + fcf + ki • k 2 + k 2 ,) 2 


d 3 k 2 
(2tt) 3 


4S(V3fci/2)S(V3fc 2 /2) 
(fc 2 + fc 2 + ki • k 2 + k 2 ) 2 


(F18) 


d 3 k 2 


2B(V3|k 2 + k!/2|/2)B(V3|k 2 - ki/2|/2) 


(F19) 


(2tt) 3 


(« 2 + fc 2 + 3fc 2 /4) 2 



n/(k x ) = 3(V3/2) 6 
n z/ (k x ) = 3(73/2) 6 
n/i/Cki) = 3(73/2) 6 
niv(ki) = 3(V3/2) 6 

We shall now take the large k\ limit, or equivalently formally the kq — > limit for a fixed ki- From the asymptotic 
behavior (|F9|) we see that B{k\) 2 involves a sum of "oscillating" terms involving fc 2s ° or k^ 2s ° , and of "non-oscillating" 
terms. We shall calculate first the resulting non-oscillating contribution, then the resulting oscillating one, up to order 
included. 



3. Non-oscillating contribution up to 0(l/kf) 

We consider the small Kq limit successively for each of the four components of n(fci) in (|F15[) . 
Contribution /: Taking directly kq — > in the integral defining nj, replacing B(ki) by its asymptotic behavior (|F9 
and averaging out the oscillating terms fc^ 2 ^ gives the leading behavior 



MM ~ -^Ag 2 -^. (F20) 



Contribution II: In the integrand of (|F17[) . we use the splitting 

(fc 2 + k 2 + ki • k 2 + K 2 )- 2 = fcf 4 + [(fc 2 + fc 2 + ki • k 2 + k 2 )- 2 - fcf 4 ] . (F21) 

The first term in the right hand side gives a contribution exactly scaling as 1/kf. In the contribution of the second 
term in the right hand side, one may take the limit kq — > and replace -B 2 (-\/3fc2/2) by its asymptotic expression to 
get the subleading 1/fcj 5 contribution. Performing the change of variable k 2 = fciq in the integral and averaging out 
the oscillating terms k l 2s ° gives 

M(*i)> = ^ " ^ W^T + °(Vfci 5 ): (F22) 

with 



y?L2B 2 {V r Ak 2 /2). (F23) 

We calculate C from the exact expression (|F8I) of B: We integrate over solid angles and we use the change of variables 
-^fc 2 = kq sinha, where a varies from zero to +oo, to take advantage of the fact that 

2tt 5 / 2 

B(k sinh a) = -AL — (e s ° a - e~ s ° a ) . (F24) 

Kq smh a cosh a v 
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This leads to 



C= 12tt 3 (v / 3/2) : 



K 



+00 2 . 

da — 



(e 2s ° Q +c.c.) 



cosh( 



(F25) 



where we used the fact that AQ, = —\AT^\ 2 . The resulting integral over a may be extended over the whole real axis 
because the integrand is an even function of a; it may then be evaluated by using the general result (that we obtained 
with contour integration) 



K(9,s) 



+00 



da ■ 



2-7T sinh(s6 



cosh a + cos 8 sin 8 sinh(s7r) 
where s is a real number and 9 G]0,7r[. One simply has to take 8 — n/2, s = and s — \sq\ respectively. We get 



(F26) 



C 



24tt 4 / \/3 \ 2sinh 2 (|s |7r/2) 



cosh(|so|7r) 



(F27) 



This, together with (|D5|) . leads to the explicit expression (|F1[) for C. 

Contribution III: We directly take the limit kq — > and we replace the factors B by their asymptotic expressions 
in (|F18[) . After the change of variable k 2 = feiq, angular integration and averaging out of the oscillating terms k^ 2s ° , 
this gives 



9 4^ 5 |A^| 2 



+ OO 



dq - 



o(l/kl). 



(F28) 



In this result, we change the integration variable setting q = e Q , where a varies from — 00 to +00. The odd component 
of the integrand (involving sinha) gives a vanishing contribution. The even component of the integrand involves a 
rational fraction of cosh a to which we apply a partial fraction decomposition. Then we use (|F26[) to obtain 



n nl (ki 



4tt 5 |AT v ,| 2 3V3sinh(7r|s |/3) + sinh(27r|s |/3) 



fef 



2tt 



sinh(7r|so| 



+ o{l/k\). 



(F29) 



Contribution IV: We directly take the limit kq — > and we replace the factors B by their asymptotic expressions 
in (|F19|) . We perform the change of variable k 2 = (fci/2)q, we average out the oscillating terms kf 2s ° . The angular 
integration in spherical coordinates of axis the direction of ki may be performed using 



dv 



1 



s /2 



so/2 



/so, 



where the variable v is restricted to the interval (—1, 1). This leads to 



(rajy(ki)) = 



47r 5 |A^,| 2 36 
k\ ir 2 





/ 9 + 


1 y° 




"^o- 1 ( 






- c.c. 


J<?- 





0(1/*;?). 



(F30) 



(F31) 



Calculating this integral directly is not straightforward because of the occurrence of the absolute value \q — 1|. We 
thus split the integration domain in two intervals. For q € [0, 1] we set q = (X — 1)/(X + 1) (an increasing function 
of X, where X spans [1, +00]). For q € [1, +00] we set q = (X + 1)/(X — 1) (a decreasing function of X, where X 
here also spans [1, +00]). Then 



(n IV (ki)) = 



4^ 5 |Ag 2 9 
fcf 2tt 2 



+00 dX (1 2 -1 + x- 2 )(x - x- 1 ) 



(X 2 



l + X 



-2\2 



L X S0 - s-'X 



1 TT-sol 



o(l/kl). 



(F32) 



We then set X = e°, where a ranges from zero to +00, and we use the fact that the resulting integrand is an even 
function of a to extend the integral over the whole real axis. We integrate by parts, integrating the factor sin(a|so|), 
and we perform a partial fraction decomposition of the resulting rational fraction of cosh a. Using (IF26[) and its 
derivatives with respect to 8, we get 



. 12tt 5 |A^| 2 -6[cosh(27r|s Q |/3)-cosh(7r|so|/3)]+V3|so|[sinh(27r|so|/3) + sinh(7r|so|/3)] 
{niv(ki)) = — x 



k\ 



27r I so I sinh(7r|so| 



■o(l/kt). 
(F33) 
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Sum of the four contributions: Summing up the terms in l/k\ of the contributions nj, nji, rim and njv, we 
obtain as a global prcfactor the quantity 

V3 cosh(27r|s |/3) - cosh(7r|s |/3) 
8 | s 1 smh(7r|s |) 

Multiplying (|D2I) on both sides by sinh(|so|7r/2) and using 

2 sinha sinh b = cosh(a + b) — cosh(a — b), Va, b (F35) 

we find that S is exactly zero. As a consequence, the non-oscillating part of the momentum distribution of an infinite 
scattering length Efimov trimer behaves at large k as 

(n(k 1 )) = ^+ (l/fcf). (F36) 



4. Oscillating contribution at large ki 

In the large k\ tail of the momentum distribution, we now include oscillating terms, having oscillating factors such 
as fc^ 2 * . The calculation techniques are the same of in the previous subsection, so that we give here directly the 
result. We find that the leading oscillating terms scale as l/k\: 



n(k 1 )-(n(k 1 )) = -^|^| 2 



A 



'k x Vf 

K 



o{l/k\) 



(F37) 



where the complex amplitude A is the sum of the contributions coming from each of the four components 
(IF16IF17IF18IF19|) of the moment distribution, 



A = Ai + An + A 



ii -f j\in -f Aiv 



We successively find 
3 



Ai 
An 
Am 



87T^ 

3 
4^2 



+00 

dq 



q 2 _ V3 
q 4 + q 2 + l 16tt' 

Jo 1 



>/3sinh(47r|s Q |/3)+sinh(27r|so|/3) 



4?r 



sinh(27r|so|) 



= 77-7 / d 1 



+ 00 



f° y/Z sinh(27r|s |/3) + sinh(7r|s |/3) - iV3[cosh(2vr|so|/3) - cosh(7r|s |/3)] 



Aiv = 



2tt 2 

12 o-2s„ 
2 



q 4 + q 2 + 1 4tt 

+ °° q(i + q 2 y° fZi/^+i 2 ) n~ v 2 yo/ 2 

' 1 dv ■ 



sinh(7r|s |) 



dq 



{q 2 + 3) 2 (q 2 + 1) 



1-v 2 



0.0243657158- 0.0698680251i. 



(F38) 

(F39) 
(F40) 
{F41) 
(F42) 



We have calculated analytically all these integrals, except for (|F42I) where the angular integration gives rise to the 
integral over v and thus to a difficult hypergeometric function. We used numerical integration for (|F42[) . Finally 



A ~ 0.1022397786- 0.1218775240i. 



(F43) 



5. Momentum distribution at the origin 

The contribution rij (ki ) is straightforward to calculate at all ki : 

ni[kl) 4™ ^ 2 ) (k 2 + 4^/3) I/ 2 ' ^ 44j 

The contribution nrj(ki) is also exactly calculable by performing the change of variable = (2/\/3) sinha and using 
the generalization of (IF26[) : 

+oo e isa 27TSill[ a (i7r-ao)] 

da „„„u I Z^T- I = ' ( F45 ) 



oo 



cosh a — cosh ao sinh ao sinh(s7r) 
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where ao is a complex number with non-zero imaginary part. This allows to obtain an exact expression of ni n (k-i) 
if one further applies integration by part, integrating the factor sin(|so|a). 

We do not give here the resulting expressions. Contrarily to these first three contributions to n(ki), the contribution 
riivQti) in (|F19I) indeed seems difficult to calculate analytically for an arbitrary kj., and blocked our attempt to 
calculate n(k x ) exactly. For ki = however it becomes equal to the contribution rin and may be evaluated exactly. 
We have thus calculated the value of n(ki = 0): 



ni(0) 



rajj(O) 



m v (0) 
with B(0) = 




3B 2 (0) 
8ttkq 

6V3B 2 (0) 
ir\s \ 2 K 

'so | cosh(|s |27r/3) 



(F46) 



+ 



cosh(|so|7r) 
cosh(|s |47r/3) 



sinh(|so|7r) cosh(|so|7r) 



5^/3 
18 



sinh(|s |47r/3) + sinh(|s |27r/3) 
sinh(|so|7r) cosh(|so|7r) 



- 2 



uizi(O) = 



6V3B 2 (0) 



7t|so|ko 
= n/z(0), 

^47T 5 /2| So | 




cosh(|s |7r/3) - cosh(|s |27r/3) + (|s |/%/3)[sinh(|s |27r/3) + sinh(|s |7r/3)] 



sinh(|s |7r) 



(F47) 

(F48) 
(F49) 



'«§ according to jF8|. This leads to (fl63)) . 



Appendix G: A direct calculation of E(A) 



To calculate the cut-off dependent energy E(A) defined in (|164j) for an infinite scattering length Efimov trimcr, 
the method consisting in calculating the momentum distribution n(k) and then integrating (|164j) is numerically 
demanding: a double integral has to be performed to obtain n(k), see (|F19|) . so that the evaluation of E(A) results 
in a triple integral. A more direct formulation, involving only a double integration, is proposed here. One simply 
rewrites (|164|) as 



E(A) 



rt z k h 2 k 2 



(27T) 



n{k) - 



(Gl) 



where the function f(k) is equal to unity for < k < A and is equal to zero otherwise. Then one plugs in (|G1|) the 
expression (|F15|) of n(k), also replacing C with its integral expression (|F23[) . An integration over two vectors in R 3 
appears: 

(G2) 
(G3) 
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easy (A J — 6 I 
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k ■ q -r n, 
k/2|) 
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{q 2 + l k 2 + k 2 q) 2 



The first part -E eaS y of this expression originates from the bits n/, n//, nm of the momentum distribution and from 
C; angular integrations may be performed, one is left with a double integral over the moduli k and q. Taking kq as a 
unit of momentum and h Kg/m as a unit of energy in what follows: 
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B 2 (^k) + 2B 2 {^q) + 4B(^k)B(^q) 2B 2 (^q) 



(k 2 + q 2 + l) 2 - k 2 q 2 



that we integrate numerically. The second part i?hard(A) in (|G4j) originates from the bit njv of the momentum 
distribution. Performing the change of variables q = (ki — k2)/2 and k = ki + k2 ensures that the factors B are now 
functions of the moduli k\ and k 2 only, 
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so that angular integrations may again be performed, involving the integral 
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(G7) 
(G8) 



where u is the cosine of the angle between the vectors and k 2 , U = [A 2 — (k\ + k 2 )]/ (2k 1 k 2 ), max(a, b) (resp. 
min(a, 6)) is the largest (resp. smallest) of the two numbers a and b, and the notation [F(u)}^ stands for Fib) — F{a) 
for any function F(u). We also used the fact that |ki + k 2 | < A if and only if u < U. This leads to 



£ hard (A) = 3 (^pj (j^) dh k\ £°° dk 2 k\ I^M^Cy^bCy^). 



(G9) 



Further simplifications may be performed. One can map the integration to the domain k\ > k 2 since the integrand is 
a symmetric function of fci and k 2 . Then performing the change of variable k\ = q + k/2 and k 2 = q — k/2, and using 
the fact that I{k\,k 2 ) = if k\ — k 2 > A, we obtain the useful form 
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that we integrate numerically. A useful result to control the numerical error due to the truncation of the integral over 
q to a value 3> A and S> 1 is 
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Appendix H: Validity of the energy-momentum relation for 77 — >■ + with a smooth regularisation 

Here we prove (|176p for the skeptics. We take hq as a unit of wavevector and h 2 K^/m as a unit of energy, so that 
the energy of the infinite scattering length bosonic Efimov trimer is -Etrim = — 1. 

First we obtain an integral expression for E v for a non-zero rj, using the same technique as in Appendix [UJ after 
having replaced B by B 1t in (|F15[) . After angular integration we obtain 
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We collect all the squared terms in B 2 , transforming B 2 (^^q^j into B 2 (j^-kj by an exchange of the integration 
variables q and k. As a consequence the integral over q can be performed explicitly for these terms: 
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Using the same exchange trick we that some simplification occurs among the crossed terms in B 2 (^^q^j B 2 (^-^kj. 
For convenience we split the final result in three pieces: 
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In a second step we use the property 
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(H7) 

(H8) 



| so | sinh a cosh a 

Hence we perform the change of variable k = (2/ y/3) sinh a and q = (2 /\/3) sinh /3. The first piece is transformed into 
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Taking the limit 77 — > + in we see that 
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The second piece is transformed into 



= 3 
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Calculation of this double integral looks hopeless. However with the natural change of variables 

y + x 



a = 



of Jacobian equal to 1/2, one can use the magic identity 
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Using the fact that the integrand is a symmetric function of a and j3, we can restrict the integration domain to a > (3 
that is x > so that, with the well known relation sin a sin b = [cos(a — b) — cos(a + b)}/2, we obtain 
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Since this integrand is now a sum of factorized terms, one of the integrals may be calculated (in some cases, one needs 
to exchange the order of integration over x and y). We are left with a single integration, in which we may take the 
limit rj -> 0+: 
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The third piece is transformed into 
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Using the change of variables (|H12[) . (|H13j) and the magic relation ()H14|) . we obtain the simpler form 
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We can take directly the limit 7/ — > + without producing diverging terms in E^ . The integrand is a sum of factorized 
terms; when a term involves the factor cos(|so|a;), we calculate the integral over y; when a term involves the factor 
cos(|so|y), we calculate the integral over x. We are thus left with single integration: 
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Finally, it remains to collect all the three pieces in E„. The terms proportional to \jr\ in £^ and i?^ 2 " 1 can be 
checked to cancel exactly: one uses integration by part to eliminate the logarithmic function in the integrand of the 
coefficient of 1/r) in E^\ and the resulting integrals in E^ and E^ may be calculated using (|F26[) . What remains 
is 
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(H20) 



This may be calculated using again (|F26[) . From the value of B(0) given in the Appendix [Fj we finally obtain the 
expected result 

lim £L = -1. (H21) 

7,->0 + 



Appendix I: Momentum distribution asymptotics of an Efimov trimer in presence of a smooth regularisation 



To understand the deviation between the true Efimov trimer energy -Btrim = — h 2 n 2 /m and the value E predicted 
by a at first sight convincing application of an energy-momentum relation, see (|164p . we suggested in the main text 
to apply the regularisation procedure (|174[) depending on a parameter rj that one eventually sets to + . Here we give 
the expressions of the coefficients C n and D v of the asymptotics (|175p of the corresponding momentum distribution 
n v (k). 

The calculations are similar to the one of Appendix [F] We shall need simply the asymptotic behavior of n n (k) after 
having averaged out the 0(1/ k 5 ) contributions involving oscillating terms in k ±2% ' s °': 
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The coefficient D n vanishes for rj —> + but it is non zero for rj > 0: 
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The contributions I v , J v and K n originate respectively from the bits nn, nm and njv in (|F15I) . Their expressions 
allow a numerical calculation of D v if desired. At first sight, the calculation of K v is more difficult because it involves 
a double integration; since the inner integral is from to a function of q, it may however be advanced step by step 
in parallel with the evaluation of the outer integral, so that the complexity remains the same as for a single integral. 
Anyway, such a numerical calculation for a finite r\ is not necessary, what matters is the knowledge of the derivative 
dD v /dr] in t\ = 0, see (|180|) . We obtain for the derivatives: 
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Remarkably, in (|I9|) a single integral is obtained, after integration by part, taking the derivative of the bit 
Jo 



dv . . .. All the integrals may be calculated numerically to a high precision with Maple, resulting in 

dD n , 
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3.3720476291291. 
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[132] The corresponding state vector is *) = [N\/(N t \Ni\)] 1/2 A ( f, ■ • • , t, 4, ■ • ■ > 4) ® W) where there are N t spins f and 

N± spins 4_i an d the operator A antisymmetrizes with respect to all particles. The wavefunction ip(ri, . . . ,rjv) is then 

proportional to ((f, . . . ,t,4, ■••,41® (ri, . . . ,rjv|) |*). 
[133] For purely attractive interaction potentials such as the square-well potential, above a critical particle number, the ground 

state is a collapsed state and the zero-range regime can only be reached for certain excited states. 
[134] Our derivation is similar to the one given in the two-body case and sketched in the many-body case in Section 3 of [62| . 
[135] In second quantization, (@JJ corresponds to g$ (R + r/2, R - r/2) = ($}(R + r/2)${(R- r/2)^ 4 (R- r/2)$ t (R + r/2)). 
[136] In 2D and 3D, our result does not follow from the well-known fact that, for a finite-range interaction potential in 

continuous space, — ^— VJ^ AG^}(r = 0) equals the kinetic energy; indeed, the Laplacian does not commute with the 

zero-range limit in that case. 

[137] In the lattice model in 3D, the coupling constant go is always negative in the zero-range limit |a| 3> b, and is an increasing 
function of — 1/a, as can be seen from (|13[) . 
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[138] Eq. ()48[) shows that dEo/d(\na) is proportional to the probability of finding two particles very close to each other, and it 
is natural that this probability decreases when one goes from the BEC limit (In a — > — oo) to the BCS limit (In a — > +00), 
i.e. when the interactions become less attractive [in the lattice model in 2D, the coupling constant go is always negative 
in the zero-range limit a> b, and is an increasing function of In a, as can be seen from (|14[) ]. 

[139] In 2D we assume, to facilitate the derivation, that V(r) — for r > b, but the result is expected to hold for any V(r) 
which vanishes quickly enough at infinity. 

[140] We consider two particles of opposite spin in a cubic box of side L with periodic boundary conditions, and we work in 
the limit where L is much larger than |a| and b. In this limit, there exists a "weakly interacting" eigenstate ip whose 
energy is given by the "mean- field" shift E = g/L 3 . The Hellmann-Feynman theorem gives go dE/dgo = Ei n t[ip]. But the 
wavefunction ^)(n,r2) — $(ri2)/L 3 where $ is the zero-energy scattering state normalized by <E> — > 1 at infinity. Thus 
Bint = Jd 3 r V(r)\$>(r)\ 2 /L 3 . The desired Eq.([T20)l then follows, since <E> = -a<j>. 

[141] In second quantization, ([15011 corresponds to g {2) (R + r/2, R— r/2) = (ft (R+ r/2)ft (R - r/2)$(R - r/2)V»(R + r/2)) . 

[142] See Chapter 3 in 

[143] The value of this constant is irrelevant for what follows. It could be calculated e.g. by equating the ener gies o f the weakly 
bound Efimov trimers of the lattice model with the ones of the zero-range model. This was done in [75l . Il08| . not for the 
lattice model, but for a Gaussian separable potential model. 

[144] The zero-range limit for a fixed Rt can be taken by repeatedly dividing b by the discrete scaling factor exp(7r/|so|) where 
so is defined in (|162|l and by adjusting go so that a remains fixed. In this limit the ground state energy tends to —00 as 
follows from the Thomas effect, but the restriction of the spectrum to any fixed energy window converges (see e.g. [75|]). 

[145] We note that it was already speculated in [63l ] that, in presence of the Efimov effect, "a three-body analog of the contact" 
may "play an important role". 

[146] The following tricks are used in the numerical calculation of the integral appearing in (|170[> . The integral is split in 
j-fcmax _|_ J+°° _ A linear scale is used to discretized k/no in between and 10. For k in between 10 and fc max (fc max 
is either 1000 or 5500), ln(fc/Ko) is discretized on a linear scale. The integral from fc max to infinity is estimated with 
the formula 8n(k) ~ — 36/(fc/fco) 6 , an approximate asymptotic expression that was tested numerically over the range 
100 < fc/fto < 1000. A simple test of this overall procedure is to check the normalization of n(k). It is found that the 
numerical calculation gives the correct normalization factor within a ~ 10~ 6 relative error. We note that calculating 8n(k) 
with a one percent error for k = 5500 requires a calculation of n(k) with a relative error ~ 2 x 10~ 10 . 

[147] The argument of the exponential in (|174p is a regular function of k around k — 0, contrarily to the more primitive choice 
e~ 2r > ln(fc/Ko) , and is quite natural considering flF8j) and (|F24|) . 

[148] This relation also breaks down for a Fermi-Fermi mixture with a mass ratio equal to the critical value 13.6 . . . above which 
the Efimov effect occurs. Indeed, the momentum distribution then has a subleading contribution 5n a (k) oc 1/fc 5 [l3lj . 
leading to a divergent integral in this relation. For a mass ratio strictly smaller than the critical ratio the integral converges, 
because 8n a (k) oc l/fc 5+2s where s > is the scaling exponent of the three-body wavefunction, f/i(Ari, Ar-2, Ar3) oc \ s ~ 2 
for A ->■ |S3. ' 

[149] For Efimovian eigenstates, computing the derivative of the energy with respect to the effective range would require to use 
a regularisation procedure similar to the one employed in free space in [87l. |89(|. However the derivative with respect to 
1/a can be computed fnfy . 

[150] Here we used the value of the effective range r e — 1.435 ro [124| for the Gaussian interaction potential V(r) = — Vbe _r 

with Vo equal to the value where the first two-body bound state appears. 
[151] For a finite-range potential one has g^2(r) = n/2— r 2 mi?km/(3ft 2 V) + . . . where V is the volume; the kinetic energy diverges 

in the zero-range limit as E kin ~ -E in t, thus E kin ~ -C"/(47r) 2 / d 3 r V [r)\<f>(r)\ 2 from so that £ kin ~ Civti 2 /(32mfc) 

for the square- well interaction. This behavior of g^ir) only holds at very short distance r <S b and is below the resolution 

of the Monte-Carlo data. 

[152] Strictly speaking, (|136[) was derived for a parabolic kinetic energy dispersion relation, whereas the dispersion relation in 
the Quantum Monte Carlo Hubbard lattice model deviates from a parabola at large k. The deviations however scale as 
k 4 r 2 at low k, so are expected to give rise to a higher order, 0(r 2 ) devi ation to T c . 

[153] Tan also derived Equation (|54[) independently from 68] using Ref. 125]'s Equation (17). 

[154] In other words, the Dirac distributions originating from the action of the Laplacian onto the 1/rij divergences can be 
ignored. 



